text:
const ComponentsByType = Dict{DataType, Dict{String, <:InfrastructureSystemsType}}

struct Components
    data::ComponentsByType
    time_series_storage::TimeSeriesStorage
    validation_descriptors::Vector
end

function Components(time_series_storage::TimeSeriesStorage, validation_descriptors=nothing)
    if isnothing(validation_descriptors)
        validation_descriptors = Vector()
    end

    return Components(ComponentsByType(), time_series_storage, validation_descriptors)
end

"""
    add_component!(
                   components::Components,
                   component::T;
                   kwargs...
                  ) where T <: InfrastructureSystemsType

Add a component.

Throws ArgumentError if the component's name is already stored for its concrete type.

Throws InvalidRange if any of the component's field values are outside of defined valid
range.
"""
function add_component!(
                        components::Components,
                        component::T;
                        skip_validation=false,
                       ) where T <: InfrastructureSystemsType
    if !isconcretetype(T)
        throw(ArgumentError("add_component! only accepts concrete types"))
    end

    if !haskey(components.data, T)
        components.data[T] = Dict{String, T}()
    elseif haskey(components.data[T], component.name)
        throw(ArgumentError("$(component.name) is already stored for type $T"))
    end

    if !isempty(components.validation_descriptors) && !skip_validation
        if !validate_fields(components, component)
            throw(InvalidRange("Invalid value"))
        end
    end

    if !skip_validation && !validate_struct(component)
        throw(InvalidValue("Invalid value for $(component)"))
    end

    # TODO: this check doesn't work during deserialization.
    #if has_forecasts(component)
    #    throw(ArgumentError("cannot add a component with forecasts: $component"))
    #end

    set_time_series_storage!(component, components.time_series_storage)
    components.data[T][component.name] = component
    return
end

"""
Removes all components from the system.
"""
function clear_components!(components::Components)
    for type_ in collect(keys(components.data))
        remove_components!(type_, components)
    end
end

"""
    remove_components!(
                       ::Type{T},
                       components::Components,
                      ) where T <: InfrastructureSystemsType

Remove all components of type T.

Throws ArgumentError if the type is not stored.
"""
function remove_components!(
                            ::Type{T},
                            components::Components,
                           ) where T <: InfrastructureSystemsType
    if !haskey(components.data, T)
        throw(ArgumentError("component $T is not stored"))
    end

    components_ = pop!(components.data, T)
    for component in values(components_)
        prepare_for_removal!(component)
    end

    @debug "Removed all components of type" T
    return values(components_)
end

"""
    remove_component!(
                      components::Components,
                      component::T,
                     ) where T <: InfrastructureSystemsType

Remove a component by its value.

Throws ArgumentError if the component is not stored.
"""
function remove_component!(
                           components::Components,
                           component::T,
                          ) where T <: InfrastructureSystemsType
    return _remove_component!(T, components, get_name(component))
end

"""
    remove_component!(
                      ::Type{T},
                      components::Components,
                      name::AbstractString,
                      ) where T <: InfrastructureSystemsType

Remove a component by its name.

Throws ArgumentError if the component is not stored.
"""
function remove_component!(
                           ::Type{T},
                           components::Components,
                           name::AbstractString,
                          ) where T <: InfrastructureSystemsType
    return _remove_component!(T, components, name)
end

function _remove_component!(
                            ::Type{T},
                            components::Components,
                            name::AbstractString,
                           ) where T <: InfrastructureSystemsType
    if !haskey(components.data, T)
        throw(ArgumentError("component $T is not stored"))
    end

    if !haskey(components.data[T], name)
        throw(ArgumentError("component $T name=$name is not stored"))
    end

    component = pop!(components.data[T], name)
    prepare_for_removal!(component)
    @debug "Removed component" T name
    return component
end

"""
    get_component(
                  ::Type{T},
                  components::Components,
                  name::AbstractString
                 )::Union{T, Nothing} where T <: InfrastructureSystemsType

Get the component of concrete type T with name. Returns nothing if no component matches.

See [`get_components_by_name`](@ref) if the concrete type is unknown.

Throws ArgumentError if T is not a concrete type.
"""
function get_component(
                       ::Type{T},
                       components::Components,
                       name::AbstractString
                      )::Union{T, Nothing} where T <: InfrastructureSystemsType
    if !isconcretetype(T)
        throw(ArgumentError("get_component only supports concrete types: $T"))
    end

    if !haskey(components.data, T)
        @debug "components of type $T are not stored"
        return nothing
    end

    return get(components.data[T], name, nothing)
end

"""
    get_components_by_name(
                           ::Type{T},
                           components::Components,
                           name::AbstractString
                          )::Vector{T} where T <: InfrastructureSystemsType

Get the components of abstract type T with name. Note that
InfrastructureSystems enforces unique names on each concrete type but not
across concrete types.

See [`get_component`](@ref) if the concrete type is known.

Throws ArgumentError if T is not an abstract type.
"""
function get_components_by_name(
                                ::Type{T},
                                components::Components,
                                name::AbstractString
                               )::Vector{T} where T <: InfrastructureSystemsType
    if !isabstracttype(T)
        throw(ArgumentError("get_components_by_name only supports abstract types: $T"))
    end

    components_ = Vector{T}()
    for subtype in get_all_concrete_subtypes(T)
        component = get_component(subtype, components, name)
        if !isnothing(component)
            push!(components_, component)
        end
    end

    return components_
end

"""
    get_components(
                   ::Type{T},
                   components::Components,
                  )::FlattenIteratorWrapper{T} where T <: InfrastructureSystemsType

Returns an iterator of components. T can be concrete or abstract.
Call collect on the result if an array is desired.

See also: [`iterate_components`](@ref)
"""
function get_components(
                        ::Type{T},
                        components::Components,
                       )::FlattenIteratorWrapper{T} where T <: InfrastructureSystemsType
    if isconcretetype(T)
        components_ = get(components.data, T, nothing)
        if isnothing(components_)
            iter = FlattenIteratorWrapper(T, Vector{Base.ValueIterator}([]))
        else
            iter = FlattenIteratorWrapper(T,
                                          Vector{Base.ValueIterator}([values(components_)]))
        end
    else
        types = [x for x in get_all_concrete_subtypes(T) if haskey(components.data, x)]
        iter = FlattenIteratorWrapper(T, [values(components.data[x]) for x in types])
    end

    @assert eltype(iter) == T
    return iter
end

"""
    iterate_components(obj) where T <: InfrastructureSystemsType

Iterates over all components.

# Examples
```julia
for component in iterate_components(obj)
    @show component
end
```

See also: [`get_components`](@ref)
"""
function iterate_components(components::Components)
    Channel() do channel
        for comp_dict in values(components.data)
            for component in values(comp_dict)
                put!(channel, component)
            end
        end
    end
end

function iterate_components_with_forecasts(components::Components)
    Channel() do channel
        for comp_dict in values(components.data)
            for component in values(comp_dict)
                if has_forecasts(component)
                    put!(channel, component)
                end
            end
        end
    end
end

function JSON2.write(io::IO, components::Components)
    return JSON2.write(io, encode_for_json(components))
end

function JSON2.write(components::Components)
    return JSON2.write(encode_for_json(components))
end

function encode_for_json(components::Components)
    # Convert each name-to-value component dictionary to arrays.
    new_components = Dict{String, Vector{<:InfrastructureSystemsType}}()
    for (data_type, component_dict) in components.data
        new_components[strip_module_name(data_type)] = [x for x in values(component_dict)]
    end

    return new_components
end

"""
Return an iterable of component types deserialized from JSON.
"""
function get_component_types_raw(::Type{Components}, raw::NamedTuple)
    return propertynames(raw)
end

"""
Return an iterable of components as NamedTuples deserialized from JSON.
"""
function get_components_raw(
                            ::Type{Components},
                            ::Type{T},
                            raw::NamedTuple,
                           ) where T <: InfrastructureSystemsType
    return getproperty(raw, Symbol(strip_module_name(string(T))))
end

function get_num_components(components::Components)
    count = 0
    for components in values(components.data)
        count += length(components)
    end
    return count
end

function clear_forecasts!(components::Components)
    for component in iterate_components_with_forecasts(components)
        clear_forecasts!(component)
    end
end

function get_forecast_initial_times(components::Components)::Vector{Dates.DateTime}
    initial_times = Set{Dates.DateTime}()
    for component in iterate_components_with_forecasts(components)
        get_forecast_initial_times!(initial_times, component)
    end

    return sort!(Vector{Dates.DateTime}(collect(initial_times)))
end

function get_forecasts_initial_time(components::Components)
    initial_times = get_forecast_initial_times(components)
    if isempty(initial_times)
        throw(ArgumentError("no forecasts are stored"))
    end

    return initial_times[1]
end

function get_forecasts_last_initial_time(components::Components)
    initial_times = get_forecast_initial_times(components)
    if isempty(initial_times)
        throw(ArgumentError("no forecasts are stored"))
    end

    return initial_times[end]
end

"""
    check_forecast_consistency(components::Components)

Throws DataFormatError if forecasts have inconsistent parameters.
"""
function check_forecast_consistency(components::Components)
    if !validate_forecast_consistency(components)
        throw(DataFormatError("forecasts have inconsistent parameters"))
    end
end

function validate_forecast_consistency(components::Components)
    # All component initial times must be identical.
    # We verify resolution and horizon at forecast addition.
    initial_times = nothing
    for component in iterate_components_with_forecasts(components)
        if !validate_forecast_consistency(component)
            return false
        end
        component_initial_times = Set{Dates.DateTime}()
        get_forecast_initial_times!(component_initial_times, component)
        if isnothing(initial_times)
            initial_times = component_initial_times
        elseif initial_times != component_initial_times
            @error "initial times don't match" initial_times, component_initial_times
            return false
        end
    end

    return true
end


text:
#
# include("../src/Mozi.jl")
# using Test
# using Logging
#
# using .Mozi
#
# const PATH=pwd()
# macro showbanner(word,total=99)
#     n=length(word)
#     m=(total-n)÷2
#     for i in 1:m
#         print("-")
#     end
#     print(word)
#     for i in 1:total-m-n
#         print("-")
#     end
#     println()
# end

@showbanner "Basic quad test"
st=Structure()

add_uniaxial_metal!(st,"steel",2e11,0.3,7849.0474)

add_node!(st,1,1,1,0)
add_node!(st,2,-1,1,0)
add_node!(st,3,-1,-1,0)
add_node!(st,4,1,-1,0)

add_quad!(st,1,4,1,2,3,"steel",1e-3)

set_nodal_restraint!(st,1,true,true,true,true,true,true)
set_nodal_restraint!(st,2,true,true,true,true,true,true)

lcset=LoadCaseSet()
add_static_case!(lcset,"MEMB",0)
add_static_case!(lcset,"BEND",0)
add_nodal_force!(lcset,"MEMB",4,0,-1e5,0,0,0,0)
add_nodal_force!(lcset,"BEND",4,0,0,-10,0,0,0)

assembly=assemble!(st,lcset,path=PATH)
solve(assembly)

r=result_nodal_displacement(assembly,"MEMB",4)
@show r
r=result_nodal_displacement(assembly,"BEND",4)
@show r
@test r≈[0.0, 0.0, -1.06318, 0.799041, 0.226482, 0.0] atol=1e-3

@showbanner "Isoparam quad test"
st=Structure()
add_uniaxial_metal!(st,"steel",2e11,0.3,7849.0474)

add_node!(st,1,7,3.5,0)
add_node!(st,2,5,3,0)
add_node!(st,3,5,1.5,0)
add_node!(st,4,6,2,0)
add_quad!(st,1,4,3,2,1,"steel",1e-3)

set_nodal_restraint!(st,1,true,true,true,true,true,true)
set_nodal_restraint!(st,2,true,true,true,true,true,true)
# for i in [3,4]
#     set_nodal_restraint!(st,i,false,false,false,false,false,true)
# end

lcset=LoadCaseSet()
add_static_case!(lcset,"MEMB",0)
add_static_case!(lcset,"BEND",0)
add_nodal_force!(lcset,"MEMB",4,0,-1e5,0,0,0,0)
add_nodal_force!(lcset,"BEND",4,0,0,-10,0,0,0)

assembly=assemble!(st,lcset,path=PATH)
solve(assembly)

r=result_nodal_displacement(assembly,"MEMB",4)
@show r
r=result_nodal_displacement(assembly,"BEND",4)
@test r≈[0.0, 0.0, -0.348873, 0.233849, 0.0296653, 0.0] atol=1e-3
#
# @showbanner "Quad cantilever test"
# st=Structure()
# lcset=LoadCaseSet()
# #
# add_uniaxial_metal!(st,"steel",2e11,0.3,7849.0474)
# add_node!(st,1,0,0,0)
# add_node!(st,2,1,0,0)
# add_node!(st,3,2,0,0)
# add_node!(st,4,3,0,0)
# add_node!(st,5,0,1,0)
# add_node!(st,6,1,1,0)
# add_node!(st,7,2,1,0)
# add_node!(st,8,3,1,0)
#
# add_quad!(st,1,1,2,6,5,"steel",1e-3)
# add_quad!(st,2,2,3,7,6,"steel",1e-3)
# add_quad!(st,3,3,4,8,7,"steel",1e-3)
# add_static_case!(lcset,"DL",0)
# add_nodal_force!(lcset,"DL",8,0,-1e5,0,0,0,0)
#
# set_nodal_restraint!(st,1,true,true,true,true,true,true)
# set_nodal_restraint!(st,5,true,true,true,true,true,true)
#
# assembly=assemble!(st,lcset,path=PATH)
# solve(assembly)
#
# r=result_nodal_displacement(assembly,"DL",8)
# @show r

ext: jl


text:
"""
Embedds all produced media in one big html file
"""
function generate_preview(media_root, path=joinpath(@__DIR__, "preview.html"))
    open(path, "w") do io
        for folder in readdir(media_root)
            media = joinpath(media_root, folder, "media")
            if !isfile(media) && ispath(media)
                medias = joinpath.(media, readdir(media))
                println(io, "
 $folder
")
                embed_media(io, medias)
            end
        end
    end
end

function tourl(path)
    if Sys.iswindows()
        # THere might be a nicer way?
        # Anyways, this seems to be needed on windows
        if !startswith(path, "http")
            path = "file:///" * replace(path, "\\" => "/")
        end
    end
    return repr(path)
end

# NOTE: `save_media` is the function you want to overload
# if you want to create a Gallery with custom types.
# Simply overloading the function should do the trick
# and ReferenceTests will take care of the rest.

function save_media(entry, x::Scene, path::String)
    path = joinpath(path, "image.png")
    save(FileIO.File(DataFormat{:PNG}, path), x) # work around FileIO bug for now
    [path]
end

function save_media(entry, x::String, path::String)
    out = joinpath(path, basename(x))
    if out != x
        mv(x, out, force = true)
    end
    [out]
end

function save_media(entry, x::AbstractPlotting.Stepper, path::String)
    # return a list of all file names
    images = filter(x-> endswith(x, ".png"), readdir(x.folder))
    return map(images) do img
        p = joinpath(x.folder, img)
        out = joinpath(path, basename(p))
        mv(p, out, force = true)
        out
    end
end

function save_media(entry, results::AbstractVector, path::String)
    paths = String[]
    for (i, res) in enumerate(results)
        # Only save supported results
        if res isa Union{Scene, String}
            img = joinpath(path, "image$i.png")
            save(FileIO.File(DataFormat{:PNG}, img), res) # work around FileIO
            push!(paths, img)
        end
    end
    paths
end

function save_media(example, events::RecordEvents, path::String)
    # the path is fixed at record time to be stored relative to the example
    epath = event_path(example, "")
    isfile(epath) || error("Can't find events for example: $(example.unique_name). Please run `record_example_events()`")
    # the current path of RecordEvents is where we now actually want to store the video
    video_path = joinpath(path, "video.mp4")
    record(events.scene, video_path) do io
        replay_events(events.scene, epath) do
            recordframe!(io)
        end
    end
    return [video_path]
end

"""
    embed_image(path::AbstractString)
Returns the html to embed an image
"""
function embed_image(path::AbstractString, alt = "")
    if splitext(path)[2] == "pdf"
        return """

        """
    end
    """

    """
end

"""
    embed_video(path::AbstractString)

Generates a html formatted string for embedding video into Documenter Markdown files
(since `Documenter.jl` doesn't support directly embedding mp4's using ![]() syntax).
"""
function embed_video(path::AbstractString, alt = "")
    """
    <video controls autoplay loop muted>
      <source src=$(tourl(path)) type="video/mp4">
      Your browser does not support mp4. Please use a modern browser like Chrome or Firefox.
    </video>
    """
end

"""
Embeds the most common media types as html
"""
function embed_media(path::String, alt = "")
    file, ext = splitext(path)
    if ext in (".png", ".jpg", ".jpeg", ".JPEG", ".JPG", ".gif", ".pdf", ".svg")
        return embed_image(path, alt)
    elseif ext == ".mp4"
        return embed_video(path, alt)
    else
        error("Unknown media extension: $ext with path: $path")
    end
end

"""
Embeds a vector of media files as HTML
"""
function embed_media(io::IO, paths::AbstractVector{<: AbstractString}, caption = "")
    for (i, path) in enumerate(paths)
        occursin("thumb", path) && continue
        println(io, """

                $(embed_media(path, caption))

        """)
    end
end


text:
# Copyright (c) 2018: Matthew Wilhelm & Matthew Stuber.
# This code is licensed under MIT license (see LICENSE.md for full details)
#############################################################################
# McCormick.jl
# A McCormick operator library in Julia
# See https://github.com/PSORLab/McCormick.jl
#############################################################################
# src/implicit_routines/implicit.jl
# Subroutines for generating relaxations of implicit functions.
#############################################################################

"""
$(SIGNATURES)

An operator that cuts the `x` object using the `y bounds` in a differentiable
or nonsmooth fashion to achieve a composite relaxation within `y`.
"""
function final_cut(x::MC{N,NS}, y::MC{N,NS}) where N
    Intv = x.Intv ∩ y.Intv
    if x.cc < y.cc
        cc = x.cc
        cc_grad::SVector{N,Float64} = x.cc_grad
    else
        cc = y.cc
        cc_grad = y.cc_grad
    end
    if x.cv > y.cv
        cv = x.cv
        cv_grad::SVector{N,Float64} = x.cv_grad
    else
        cv = y.cv
        cv_grad = y.cv_grad
    end
    x_out::MC{N,NS} = MC{N,NS}(cv, cc, Intv, cv_grad, cc_grad, x.cnst && y.cnst)
    return x_out
end
function final_cut(x::MC{N,Diff}, y::MC{N,Diff}) where N
    x_out = intersect(x, y)
    return x_out
end

"""
$(TYPEDEF)

An abstract type for each manner of contractor using in the implicit function relaxation algorithms.
"""
abstract type AbstractContractorMC end

"""
$(TYPEDEF)

The Gauss-Seidel implementation of the Newton contractor used in the implicit relaxation scheme.
"""
struct NewtonGS <: AbstractContractorMC end

"""
$(TYPEDEF)

The componentwise implementation of the Krawczyk contractor used in the implicit relaxation scheme.
"""
struct KrawczykCW <: AbstractContractorMC end

"""
$(TYPEDEF)

An abstract type for each manner of preconditioner used in the implicit function relaxation algorithms.
"""
abstract type AbstractPreconditionerMC end

"""
$(SIGNATURES)

Creates storage corresponding to `x::AbstractPreconditionerMC` and `t::T where T<:RelaxTag`.
"""
function preconditioner_storage(x::AbstractPreconditionerMC, t::T) where T <: RelaxTag
    error("Must define function that generates appropriate storage type for preconditioner")
end

"""
$(TYPEDEF)

An abstract type for each manner of callback functions used in the implicit function relaxation algorithms.
"""
abstract type AbstractMCCallback end

"""
$(TYPEDEF)

A structure used to compute implicit relaxations.

$(TYPEDFIELDS)
"""
struct MCCallback{FH <: Function, FJ <: Function, C <: AbstractContractorMC,
                  PRE <: AbstractPreconditionerMC, N, T <: RelaxTag,
                  AMAT <: AbstractMatrix} <:  AbstractMCCallback
    "Function h(x,p) = 0 defined in place by h!(out,x,p)"
    h!::FH
    "Jacobian of h(x,p) w.r.t x"
    hj!::FJ
    ""
    H::Vector{MC{N,T}}
    ""
    J::AMAT
    xmid::Vector{Float64}
    X::Vector{Interval{Float64}}
    P::Vector{Interval{Float64}}
    "State space dimension"
    nx::Int
    "Decision space dimension"
    np::Int
    "Convex combination parameter"
    lambda::Float64
    "Tolerance for interval equality"
    eps::Float64
    "Number of contractor steps to take"
    kmax::Int
    pref_mc::Vector{MC{N,T}}
    p_mc::Vector{MC{N,T}}
    x0_mc::Vector{MC{N,T}}
    x_mc::Vector{MC{N,T}}
    xa_mc::Vector{MC{N,T}}
    xA_mc::Vector{MC{N,T}}
    aff_mc::Vector{MC{N,T}}
    z_mc::Vector{MC{N,T}}
    "Type of contractor used in implicit relaxation routine."
    contractor::C
    "Preconditioner used in the implicit relaxation routine."
    preconditioner::PRE
    "Boolean indicating that the preconditioner should be applied"
    apply_precond::Bool
    param::Vector{Vector{MC{N,T}}}
end
function MCCallback(h!::FH, hj!::FJ, nx::Int, np::Int,
                    contractor::S = NewtonGS(),
                    preconditioner::T = DenseMidInv(zeros(Float64,1,1), zeros(Interval{Float64},1), 1, 1),
                    relax_tag::TAG = NS()) where {FH <: Function,
                                                  FJ <: Function,
                                                  S <: AbstractContractorMC,
                                                  T,
                                                  TAG <: RelaxTag}
    H = zeros(MC{np,TAG}, (nx,))
    xmid = zeros(Float64, (nx,))
    P = zeros(Interval{Float64}, (np,))
    X = zeros(Interval{Float64}, (nx,))
    lambda = 0.5
    eps = 0.0
    kmax = 2
    p_ref = zeros(MC{np,TAG}, (np,))
    p_mc = zeros(MC{np,TAG}, (np,))
    x0_mc = zeros(MC{np,TAG}, (nx,))
    x_mc = zeros(MC{np,TAG}, (nx,))
    xa_mc = zeros(MC{np,TAG}, (nx,))
    xA_mc = zeros(MC{np,TAG}, (nx,))
    aff_mc = zeros(MC{np,TAG}, (nx,))
    z_mc = zeros(MC{np,TAG}, (nx,))
    contractor = NewtonGS()
    preconditioner = preconditioner(h!, hj!, nx, np)
    J = preconditioner_storage(preconditioner, relax_tag)
    apply_precond = true
    param = fill(zeros(MC{np,TAG}, (nx, )), (kmax,))

    return MCCallback{FH, FJ, NewtonGS, DenseMidInv, np, TAG, typeof(J)}(h!, hj!, H, J, xmid, X, P, nx, np,
                                                                         lambda, eps, kmax, p_ref, p_mc, x0_mc, x_mc,
                                                                         xa_mc, xA_mc, aff_mc, z_mc,
                                                                         contractor, preconditioner,
                                                                         apply_precond, param)
end
function (d::MCCallback)()
    d.h!(d.H, d.z_mc, d.p_mc)
    d.hj!(d.J, d.aff_mc, d.p_mc)
    return
end

include("preconditioner/dense.jl")

"""
$(SIGNATURES)

Computates the affine relaxations of the state variable.
"""
function affine_exp!(x::S, p::Vector{MC{N,T}}, d::MCCallback) where {S, N, T<:RelaxTag}

    S1 = zero(MC{N,T})
    S2 = zero(MC{N,T})
    S3 = zero(MC{N,T})
    @inbounds for i = 1:d.nx
        S1 = zero(MC{N,T})
        S2 = zero(MC{N,T})
        S3 = zero(MC{N,T})
        @inbounds for j = 1:N
            S1 += (p[j]-d.pref_mc[j])*x[i].cv_grad[j]
            S2 += (p[j]-d.pref_mc[j])*x[i].cc_grad[j]
            S3 += (d.lambda*x[i].cv_grad[j]+(1.0-d.lambda)*x[i].cc_grad[j])*(p[j]-d.pref_mc[j])
        end
        temp1 = x[i].cv + S1
        temp2 = x[i].cc + S2
        temp3 = d.lambda*x[i].cv+(1.0-d.lambda)*x[i].cc+S3
        d.xa_mc[i] = MC{N,T}(temp1.cv, temp1.cc, temp1.Intv, x[i].cv_grad, x[i].cv_grad, S1.cnst)
        d.xA_mc[i] = MC{N,T}(temp2.cv, temp2.cc, temp2.Intv, x[i].cc_grad, x[i].cc_grad, S2.cnst)
        d.z_mc[i] = MC{N,T}(temp3.cv, temp3.cc, temp3.Intv,
                            d.lambda*x[i].cv_grad+(1.0-d.lambda)*x[i].cc_grad,
                            d.lambda*x[i].cv_grad+(1.0-d.lambda)*x[i].cc_grad, S3.cnst)
    end
    return
end

"""
$(SIGNATURES)

Corrects the relaxation of the state variable `x_mc` if the affine relaxation,
"""
function correct_exp!(d::MCCallback{FH,FJ,C,PRE,N,T}) where {FH <: Function,
                                                             FJ <: Function,
                                                             C, PRE, N,
                                                             T<:RelaxTag}
    zero_grad = zeros(SVector{N,Float64})
    @inbounds for i = 1:d.nx
        if (d.z_mc[i].Intv.lo - d.eps < d.X[i].lo) && (d.z_mc[i].Intv.hi + d.eps > d.X[i].hi)
            d.x_mc[i] = MC{N,T}(d.X[i])
        elseif d.z_mc[i].Intv.lo - d.eps < d.X[i].lo
            d.x_mc[i] = MC{N,T}(d.X[i].lo, d.x_mc[i].cc, Interval(d.X[i].lo, d.x_mc[i].Intv.hi),
                                zero_grad, d.x_mc[i].cc_grad, d.x_mc[i].cnst)
        else
            if d.z_mc[i].Intv.hi + d.eps > d.X[i].hi
                d.x_mc[i] = MC{N,T}(d.x_mc[i].cv, d.X[i].hi, Interval(d.x_mc[i].Intv.lo, d.X[i].hi),
                                d.x_mc[i].cv_grad, zero_grad, d.x_mc[i].cnst)
            end
        end
    end
    return
end

include("contract.jl")

"""
$(SIGNATURES)
"""
function precond_and_contract!(callback!::MCCallback{FH,FJ,C,PRE,N,T}) where {FH <: Function,
                                                                              FJ <: Function,
                                                                              C <: AbstractContractorMC,
                                                                              PRE <: AbstractPreconditionerMC,
                                                                              N, T<:RelaxTag}
    @. callback!.aff_mc = MC{N,T}(cv(callback!.xa_mc), cc(callback!.xA_mc))
    callback!()
    if callback!.apply_precond
        precondition!(callback!.preconditioner, callback!.H, callback!.J)
    end
    contract!(callback!.contractor, callback!)
    return
end

"""
$(SIGNATURES)

Populates `x_mc`, `xa_mc`, `xA_mc`, and `z_mc` with affine bounds.
"""
function populate_affine!(d::MCCallback{FH,FJ,C,PRE,N,T}, interval_bnds::Bool) where {FH <: Function,
                                                                                FJ <: Function,
                                                                                C <: AbstractContractorMC,
                                                                                PRE <: AbstractPreconditionerMC,
                                                                                N, T<:RelaxTag}
    if interval_bnds
        @inbounds for i in 1:d.nx
            d.x_mc[i] = MC{N,T}(d.X[i].lo, d.X[i].hi)
            d.xa_mc[i] = MC{N,T}(d.X[i].lo, d.X[i].lo)
            d.xA_mc[i] = MC{N,T}(d.X[i].hi, d.X[i].hi)
            d.z_mc[i] = d.lambda*d.xa_mc[i] + (1.0 - d.lambda)*d.xA_mc[i]
         end
    end
    return
end

"""
$(SIGNATURES)

Constructs parameters need to compute relaxations of `h`.
"""
function gen_expansion_params!(d::MCCallback, interval_bnds::Bool = true) where {N, T<:RelaxTag}
    populate_affine!(d, interval_bnds)
    @. d.param[1] = d.x_mc
    for k = 2:d.kmax
        precond_and_contract!(d)
        affine_exp!(d.x_mc, d.pref_mc, d)
        correct_exp!(d)
        @. d.param[k] = d.x_mc
    end
    return
end

"""
$(SIGNATURES)

Compute relaxations of `x(p)` defined by `h(x,p) = 0` where `h` is specifed as `h(out, x, p)`.
"""
function implicit_relax_h!(d::MCCallback, interval_bnds::Bool = true) where {N, T<:RelaxTag}
    populate_affine!(d, interval_bnds)
    for k = 1:d.kmax
        affine_exp!(d.param[k], d.p_mc, d)
        precond_and_contract!(d)
    end
    return
end

max_stars_repo_name: PSORLab/RSActivationFunctions

max_issues_repo_name: PSORLab/RSActivationFunctions

max_forks_repo_name: PSORLab/RSActivationFunctions


text:
using Dash
using DashHtmlComponents
using DashCoreComponents

app = dash()

app.layout = html_div() do
    dcc_input(id="input", value = 1),
    dcc_input(id="output-1"),
    dcc_input(id="output-2"),
    dcc_input(id="output-3"),
    dcc_input(id="output-4")
end

callback!(
    ClientsideFunction("clientside", "add_to_four_outputs"),
    app, callid"input.value=>output-1.value, output-2.value, output-3.value, output-4.value"
    )

run_server(app)


text:
check_singleton_type(x) = x
check_singleton_type(x::T) where T <:ForwardDiff.Dual = x.value
function check_tuple_types(X)
    return NamedTuple{keys(X)}(check_singleton_type.(values(X)))
end

max_stars_repo_name: sqwayer/AnimalBehavior

max_issues_repo_name: sqwayer/AnimalBehavior


text:
module RunTime

using Dates: now, canonicalize, CompoundPeriod

macro runtime(args...)

    println("-"^80)

    st = now()

    @show args

    println("Starting at ", st)

    eval(args...)

    en = now()

    println("Ended at ", en)

    println("Took ", canonicalize(CompoundPeriod(en - st)))

    return println()

end

end


text:
import Gtk, Clipper, GeometryTypes, Reactive, Statistics
const gtk = Gtk
const cl = Clipper
const gt = GeometryTypes
const re = Reactive
const stat = Statistics

"""
    init_text!(config::Config, key::AbstractString, init="")

Registers a field in the gui for pritning arbitrary values.
The type of `init` determines the type, that can passed to `show_text!`.
All registered texts must have a unique `key`.
"""
function init_text!(config::Config, key::AbstractString, init="")
    config.printables[key] = re.Signal(init)
end

"""
    show_text!(config::Config, key::AbstractString, value)

Shows the given `value` in the registered field in the gui (see `init_text!`).
Uses the key for the right field.
"""
function show_text!(config::Config, key::AbstractString, value)
    push!(config.printables[key], value)
end

"""
    get_text!(config::Config, key::AbstractString)

Returns the current `value` in the registered field in the gui (see `init_text!`).
Uses the key for the right field.
"""
function get_text!(config::Config, key::AbstractString)
    re.value(config.printables[key])
end

"""
    show_graphic!(config::Config, key::AbstractString, position, graphic::GraphicElement)

Shows a GraphicElement in the gui. There are various graphic element like points, circle, line and arrows or text.
Some need mor than a position, like the circle needs also a radius. Look in the specific docs of the graphic element.
The keys need to be unique for each graphic element. Graphic elements don't need to registered beforhead.
```
"""
function show_graphic!(config::Config, key::AbstractString, position, graphic::Function)
    config.drawables[key] = (position, graphic)
end

"""
    get_manip_valu(config::Config, key::AbstractString)

Returns the current value of a manipulator by its key. For example `true` or `false` if a checkbutton is checked or not.
"""
function get_manip_value(config::Config, key::AbstractString)
    get_value(config.manipulators[key][1], config.manipulators[key][2])
end

"""
    remove_graphic!(config::Config, key::AbstractString)

Removes a graphic element from the gui by its key.
"""
function remove_graphic!(config::Config, key::AbstractString)
    delete!(config.drawables, key)
end

"""
    remove_graphics_filter!(filter::Function, config::Config)

Removes graphic elements, for which key the filter function returns true.

# Examples

```jldoctest
julia> remove_graphics_filter!(key -> startswith(key, "arrows_"), config)

```
"""
function remove_graphics_filter!(filter::Function, config::Config)
    for s in keys(config.drawables)
        if filter(s)
            delete!(config.drawables, s)
        end
    end
end

"""
    is_in_line_of_sight(ws::WorldState, my_position::AbstractVector{T}, p_view::AbstractVector{T}; tol::Real=0.0, excludes::Vector{Robot{T}}=Robot[]) where T<:Number

Check if in the current WorldState `ws` a position `p_view` is visible and there is no obstacle in the line my position->`p_view`, that prevents the sight.

# Current limits

- only the point `p_view` will be checked, not the object on this position

# Optional Arguments

- `tol::Number=0.0`: tolerace of the radii for all elements in `ws`
(the tolerance will be subtracted from the real radius)
- `excludes::AbstractVector{<:AbstractVector{<:Number}}=Vector{Float64}[]`: elemtent which are subtracted from the list in `ws`

# Examples

```jldoctest
julia> ws = WorldState([[3.0, 3.0]], [[6.0, 3.0]], [8.0, 3.0])
WorldState(Array{Float64,1}[[3.0, 3.0]], Array{Float64,1}[[6.0, 3.0]], [8.0, 3.0], 0.15, 0.05)

julia> is_in_line_of_sight(ws, ws.myPlayers[1], ws.ball)
false
```
"""
function is_in_line_of_sight(ws::WorldState, my_position::AbstractVector{T}, p_view::AbstractVector{T}; tol::Real=0.0, excludes::Vector{Robot{T}}=Robot{Float64}[]) where T<:Number
    view_line = gt.LineSegment(gt.Point2(my_position), gt.Point2(p_view))
    !any([any(intersection(view_line, gt.Circle(gt.Point2(x), get_size(x) - tol), only_in_lineSegment = true)[[1,3]]) for x in vcat(ws.myPlayers, ws.opponents) if x ∉ excludes])
end

"""
    iswithin_hfov(player::Robot{<:Number}, point::AbstractVector{<:Number})

Checks if a point is within the horizontal field of view of the player. Doesn't mean, that the robot can see the point!

# Examples

```jldoctest
julia> iswithin_hfov(Robot([2.0, 2.0], 0.0), [4.0, 2.0])
true
```
"""
function iswithin_hfov(player::Robot{<:Number}, point::AbstractVector{<:Number})
    p = point .- player
    angle_diff = abs(player.body_orientation - atan(p[2], p[1]))
    angle_diff < player.hfov/2 || angle_diff > 2π - player.hfov/2
end

"""
    add_slider!(callback::Function, config::Config, key::AbstractString, range::Range{<:T}, init::T) where T<:Number

Adds a `GtkSlider` to the gui. You can access it's current value by `get_manip_value`
"""
function add_slider!(callback::Function, config::Config, key::AbstractString, range::AbstractRange{T}, init::T) where T<:Number
    speed = gtk.GtkScale(false, range)
    speed_adj = gtk.GtkAdjustment(speed)
    gtk.set_gtk_property!(speed_adj, :value, init)
    config.manipulators[key] = (speed, T, callback)
    return
end
add_slider!(config::Config, key::AbstractString, range::AbstractRange{T}, init::T) where T<:Number = add_slider!(_ -> nothing, config, key, range, init)

"""
    check_button(name::String)

Adds a `GtkCheckButton` to the gui. You can access it's current value by `get_manip_value`
"""
function add_check_button!(callback::Function, config::Config, key::AbstractString, name::AbstractString)
    config.manipulators[key] = (gtk.GtkCheckButton(name), Bool, callback)
    return
end
add_check_button!(config::Config, key::AbstractString, name::AbstractString) = add_check_button!(_ -> nothing, config, key, name)

"""
    radio_button(name::String)

Adds and returns a `GtkRadioButton` to the gui. The radiobutton will not shown in the gui, only the callback function
will be registered. For showing the radio button, you need to put the in a radio button group. See `add_radio_button_group!`.
"""
function add_radio_button!(callback::Function, config::Config, key::AbstractString, name::String)
    rb = gtk.GtkRadioButton(name)
    config.manipulators[key] = (rb, Bool, callback)
    return rb
end
add_radio_button!(config::Config, key::AbstractString, name::AbstractString) = add_radio_button!(_ -> nothing, config, key, name)

"""
    radio_button_group(elements::AbstractVector)

Adds a `GtkRadioButtonGroup` to the gui. You must pass a vector of radio buttons to it (see `add_radio_button`).
You can access it's current choosen radio button by `get_manip_value`
"""
function add_radio_button_group!(config::Config, key::AbstractString, elements::AbstractVector)
    config.manipulators[key] = (gtk.GtkRadioButtonGroup(elements, 1), String, _ -> nothing)
    return
end

"""
Easier access to values stored in this container.
"""
get_value(widget::gtk.GtkScaleLeaf, t::Type) = gtk.get_gtk_property(gtk.GtkAdjustment(widget), :value, t)
get_value(widget::gtk.GtkCheckButtonLeaf, t::Type) = gtk.get_gtk_property(widget, :active, t)
get_value(widget::gtk.GtkRadioButtonGroupLeaf, t::Type) = gtk.get_gtk_property(gtk.get_gtk_property(widget, :active), :label, t)

"""
    intervall(start::Int, stop::Int, mod::Int)

Return a array with elements between `start` and `stop` in a modulo counter.
`start` is exclusive, `stop` inclusive. Assumes the input is 1 based.
The modulocounter is starting by 1 and counts until `mod`.

# Examples

```jldoctest
julia> intervall(3, 1, 4)
2-element Array{Int64,1}:
 4
 1
"""
function intervall(start::Int, stop::Int, mod::Int, full::Bool=false)
    if full && start==stop
        return circshift(1:mod, -start)
    end

    # initialization, correct to 0 based
    intervall_array = Int[]
    start -= 1
    stop -= 1

    # modulo count
    while start != stop
        start = (start+1) % mod
        push!(intervall_array, start)
    end

    # correction from 0 based to 1 based
    return intervall_array.+1
end

"""
    get_edge(to::AbstractVector{<:Number}, tol::Real=eps())

Checks if a 2D point `to` lies on a edge of the (hardcoded) filed within a `tol`.
By default, tol is set to `eps()`.
The edges in the world coordinate systems are named like:

(0.0,6.0) (9.0,6.0)
    +----3----+
    |         |
    4         2
    |         |
    +----1----+
(0.0,0.0) (9.0,0.0)
# Examples

```jldoctest
julia> get_edge([9.0, 2.32])
2
"""
function get_edge(to::AbstractVector{<:Number}, tol::Real=eps())
    @inbounds if isapprox(to[2], 0.0, atol=tol)
        return 1
    elseif isapprox(to[1], 9.0, atol=tol)
        return 2
    elseif isapprox(to[2], 6.0, atol=tol)
        return 3
    else
        return 4
    end
end

"""
    raster(from::Real, to::Real, n::Int)

Rasters a range [`from`, `to`] in `n` segments.
[`from`, `to`] are exclusive.

# Examples

```jldoctest
julia> raster(0.0, 2.0, 3)
[0.5, 1.0, 1.5]
"""
function raster(from::Real, to::Real, n::Int)
    padding = (to - from)/(n+1)
    range(from+padding, stop=to-padding, length=n)
end

"""
    get_shadow_polygon(position::AbstractVector{<:Number}, obstacle::Robot{<:Number}, config::Config)

Calculates the shadow of a virtial omnidirectional lightsource on the `position` for a robot obstacle. Assumes that the robot is a circle with `position` and radius `size`.

# Examples

```jldoctest
julia> config = SearchBall.get_config(["--add_myPlayer","\"--position 4.0 4.0\""]);
julia> get_shadow_polygon([0.0, 0.0], config.ws.myPlayers[1], config)
5-element Array{GeometryTypes.Point{2,Float64},1}:
 [4.10322, 3.89116]
 [0.0, 0.0]
 [0.0, 0.0]
 [3.89116, 4.10322]
 [4.10322, 3.89116]
"""
function get_shadow_polygon(position::AbstractVector{<:Number}, obstacle::Robot{<:Number}, config::Config)
        # get left (tan2) and right (tan1) tangents
        tan1 = outer_tangent(gt.Point2(position), gt.Circle(gt.Point2(obstacle), get_size(obstacle)))
        isinside1, isec1 = intersection(tan1, gt.HyperRectangle(config.field_corners))
        tan2 = outer_tangent(gt.Circle(gt.Point2(obstacle), get_size(obstacle)), gt.Point2(position), true)
        isinside2, isec2 = intersection(tan2, gt.HyperRectangle(config.field_corners))

        !(isinside1 && isinside2) && return Vector{Float64}[]

        # add corners of the field to the polygones, if the corner belongs to it
        edge1 = get_edge(isec2, 1e-10)
        edge2 = get_edge(isec1, 1e-10)
        vcat([tan2[2], isec2], config.field_corners[intervall(edge1, edge2, length(config.field_corners))], [isec1, tan1[2], tan2[2]])
end

"""
    get_shadow_polygons(current_player::Robot{<:Number}, config::Config, with_hfov::Bool=false, pos::AbstractVector{<:Number}=current_player, orientation::Number=current_player.body_orientation)

Returns an array with all shadow polygons occur by a view position `new_pos` within
the current world state `ws`.

# Examples

```jldoctest
julia> config = get_config(["--add_myPlayer","\"--position 4.0 4.0\"","--add_opponent","\"--position 6.0 4.0\""])
julia> get_shadow_polygons(config.ws.myPlayers[1], config)
1-element Array{Array{GeometryTypes.Point{2,Float64},1},1}:
 GeometryTypes.Point{2,Float64}[[9.0, 4.3761], [5.9887, 4.1496], [5.9887, 3.8504], [9.0, 3.6239]]
"""
function get_shadow_polygons(current_player::Robot{<:Number}, config::Config, with_hfov::Bool=false, pos::AbstractVector{<:Number}=current_player, orientation::Number=current_player.body_orientation)
    # TODO: only obstacle which are within line of sight and in the hfov
    polygons = Vector{gt.Point2{Float64}}[get_shadow_polygon(pos, obstacle, config) for obstacle in vcat(config.ws.myPlayers, config.ws.opponents) if obstacle != current_player]

    if with_hfov
        init_orientation = [1.0, 0.0]
        d1 = rotate2d(init_orientation, current_player.body_orientation + current_player.hfov/2)
        isinside1, isec1 = intersection(gt.LineSegment(gt.Point2(pos), gt.Point2(d1 .+ pos)), gt.HyperRectangle(config.field_corners))
        d2 = rotate2d(init_orientation, current_player.body_orientation - current_player.hfov/2)
        isinside2, isec2 = intersection(gt.LineSegment(gt.Point2(pos), gt.Point2(d2.+ pos)), gt.HyperRectangle(config.field_corners))
        if isinside1 && isinside2
            edge1 = get_edge(isec1, 1e-10)
            edge2 = get_edge(isec2, 1e-10)
            corner_numbers = intervall(edge1, edge2, length(config.field_corners))
            outofview_field = vcat([pos, isec1], config.field_corners[intervall(edge1, edge2, length(config.field_corners), current_player.hfov < π)], [isec2])
            push!(polygons, outofview_field)
        end
    end

    return polygon_union(polygons)
end

using Distributed
using SharedArrays

"""
    sample(config::Config, off_idx::Int, my_reduce::Function=(x -> polygon_area(@view(x[2:end]))), my_map::Function=get_shadow)

# Optional Arguments
- `my_reduce::Function=(x -> polygon_area(@view(x[2:end])))`: function for calculation the scalar value on the current position, by default the area of the shadow (polygon)
- `my_map::Function=get_shadow`: function, which shall be evaluation on every position on the field, by default the shadow of the current player of myPlayers

Calculates for a raster on the field the shadow area.
"""
function sample(config::Config, off_idx::Int, use_hfov::Bool=false, my_map::Function=get_shadow_polygons, my_reduce::Function=polygon_area; resolution::AbstractArray{Int}=[90, 60])
    ws = config.ws
    field_corners = config.field_corners
    # initialization
    scalar_field = SharedArray(zeros(Float64, resolution...))
    # insides = falses(resolution...)
    insides = SharedArray{Bool}(resolution...)

    x_raster = raster(field_corners[1][1], field_corners[3][1], resolution[1])
    # x_length = length(x_raster)
    y_raster = raster(field_corners[1][2], field_corners[3][2], resolution[2])
    @sync @distributed for ind_x in eachindex(x_raster)
        w = x_raster[ind_x]
        @distributed for ind_y in eachindex(y_raster)
            h = y_raster[ind_y]
            shadow_area = 0.0
            # checks if no player is at this exact position
            insides[ind_x, ind_y] = !any(map(p -> is_inside([w, h], gt.Circle(gt.Point2(p), get_size(p)*1.5+eps())), [x for x in vcat(ws.myPlayers, ws.opponents) if x!=ws.myPlayers[off_idx] && is_inside(x, field_corners)]))
            if insides[ind_x, ind_y]
                polygons = my_map(ws.myPlayers[off_idx], config, use_hfov, [w, h])
                for polygon in polygons
                    shadow_area += my_reduce(polygon)
                end
            end
            scalar_field[ind_x, ind_y] = shadow_area
            # config.set_progress_bar(ind_x/x_length)
        end
    end

    scalar_field, insides
end

using ColorSchemes: inferno
const color_scheme = deepcopy(inferno)
reverse!(color_scheme.colors)

using Printf

const color_bar = get(inferno, reshape(collect(0.0:0.01:1.0), 1, 101), :clamp)

"""
    send_background(vals::AbstractArray{<:Real}, ins::AbstractArray{Bool})

Samples the shadow function on a raster of the field and draws the scalarfield in the gui.
"""
function send_background(vals::AbstractArray{<:Real}, ins::AbstractArray{Bool}, config::Config, bounds::Union{Nothing, Tuple{<:Real, <:Real}}=nothing; allways_repaint::Bool=false)

    if typeof(bounds) <: Nothing
        # clean up the positions where the player stay (players of myPlayers can't go there)
        @inbounds mean_val = stat.median(vals[ins])
        @inbounds vals[.!ins] = mean_val
        # vals[isnan.(vals)] = mean_val
        min_val = minimum(vals)
        max_val = maximum(vals)
    else
        min_val, max_val = bounds
    end

    # make a colored image from the scalar field
    image = get(color_scheme, vals, (min_val, max_val)) # :extrema

    # draw stuff
    config.drawables["_fieldBackground"] = ([[0.0, 6.0], [9.0, 0.0]], () -> (image, allways_repaint))
    config.drawables["_color_bar"] = ([[9.2, 6.0], [9.4, 0.0]], () -> (color_bar, false))
    config.drawables["color_bar_edge"] = ([[9.2, 6.0], [9.2, 0.0], [9.4, 0.0], [9.4, 6.0], [9.2, 6.0]], x -> GraphicLine(x, fill=false))
    config.drawables["text_min"] = ([[9.6, 6.0]], x -> GraphicText(x[1], @sprintf("%.2fm²", min_val), textsize=12.0, color=[0.0, 0.0, 0.0, 1.0]))
    config.drawables["text_33p"] = ([[9.6, 4.0]], x -> GraphicText(x[1], @sprintf("%.2fm²", min_val+(max_val-min_val)/3.0), textsize=12.0, color=[0.0, 0.0, 0.0, 1.0]))
    config.drawables["text_66p"] = ([[9.6, 2.0]], x -> GraphicText(x[1], @sprintf("%.2fm²", min_val+2*(max_val-min_val)/3.0), textsize=12.0, color=[0.0, 0.0, 0.0, 1.0]))
    config.drawables["text_max"] = ([[9.6, 0.0]], x -> GraphicText(x[1], @sprintf("%.2fm²", max_val), textsize=12.0, color=[0.0, 0.0, 0.0, 1.0]))
end

"""
    clear_background!(config::Config)

# Examples

```jldoctest
julia> clear_background!()

```
"""
function clear_background!(config::Config)
    delete!(config.drawables, "_fieldBackground")
    delete!(config.drawables, "_color_bar")
    delete!(config.drawables, "color_bar_edge")
    delete!(config.drawables, "text_min")
    delete!(config.drawables, "text_33p")
    delete!(config.drawables, "text_66p")
    delete!(config.drawables, "text_max")
end


text:
module NeuralNetworks

export nn_model

using Statistics
using Flux
using OceanParameterizations.DataWrangling

function cb(train,loss)
    f() = @info mean([loss(train[i][1], train[i][2]) for i in 1:length(train)])
    return f
end

# Return a trained neural network model
function nn_model(; 𝒱=nothing, model=nothing, optimizers=nothing)
    loss(x,y) = Flux.Losses.mse(model(x), y)

    for opt in optimizers
        @info opt
        Flux.train!(loss, params(model), 𝒱.training_data, opt, cb = Flux.throttle(cb(𝒱.training_data,loss), 2))
    end

    return model
end

end #module


text:
using BlackBoxOptim

p = BlackBoxOptim.example_problems["Griewank"]

BlackBoxOptim.repeated_bboptimize(2, p, 30, [
  :generating_set_search,
  #:random_search,
  #:separable_nes,
  :adaptive_de_rand_1_bin_radiuslimited],
  30.0, 1e-6, {:SaveFitnessTraceToCsv => true})