text:
"
SelectionOperator represents the method that selects a portion of the population
for crossover and mutation.

# Fields
- `_method::Function`: selection method between individuals.
- `_nSelectedParents::Real`: number of selected individuals. Its values can either
    be an integer number, which directly represents the number of individuals that
    are going to be picked for selection, or a floating point number, which then
    represents the percentage of the population that is going to be selected.
- `needsComparison::Bool = true`: tells wether the user wants the system to compare
    between the individuals in order to select between them.
- `_samplingWithRep::Bool`: tells wether the user wants the system to perform the
    selection with replacement (`true`) or not (`false`).
- `_varArgs::Array{Any}`: arguments of the selection method, if necessary.
"
struct SelectionOperator
    _method::Function
    _nSelectedParents::Real
    _needsComparison::Bool
    _samplingWithRep::Bool
    _varArgs::Array{Any}
end # struct

"""
    getMethod(selector::SelectionOperator)::Function

Returns the method for select the parents.
"""
getMethod(selector::SelectionOperator)::Function = selector._method
# function

"""
    needsComparison(selector::SelectionOperator)::Bool

Indicates wether the selection method needs to compare between individuals
or not.
"""
needsComparison(selector::SelectionOperator)::Bool = selector._needsComparison
# function

"""
    getDefaultSelector()::Tuple{Function, Array}

Returns the default selection method, which is [`tournamentSelector`](@ref) between
4 individuals.
"""
function getDefaultSelector()::Tuple{Function, Array}
    return tournamentSelector, [4]
end # function

"""

    getNParents(selector::SelectionOperator)::Real

Returns the attribute that represents the number of parents that are going to be selected.
"""
getNSelectedParents(selector::SelectionOperator)::Real = selector._nSelectedParents
# function

"""

    samplingWithReplacement(selector::SelectionOperator)::Bool

"""
samplingWithReplacement(selector::SelectionOperator)::Bool = selector._samplingWithRep
# function

"""

    getFunctionArgs(selector::SelectionOperator)::Array

Obtains the aditional arguments associated to selector method.
"""
getFunctionArgs(selector::SelectionOperator)::Array = selector._varArgs
# function

"""
    selectParentsWithReplacement(selector::SelectionOperator, population::Array{Individual},
                                evaluator::Evaluator, rng::Random.AbstractRNG)::Array{Individual}

Performs a selection between individuals with replacement using the comparison
function of the `Evaluator`.
"""
function selectParentsWithReplacement(selector::SelectionOperator, population::Array{Individual},
                                      evaluator::Evaluator, rng::Random.AbstractRNG)::Array{Individual}

    nSelectedParents = getNSelectedParents(selector)
    compare = getCompareFunction(evaluator)
    popSize = length(population)

    parents = Array{Individual}(undef, nSelectedParents)

    method = getMethod(selector)

    functionArgs = getFunctionArgs(selector)

    for i=eachindex(parents)
        if needsComparison(selector)
            index = method(population, compare, rng, functionArgs...)
        else
            index = method(population, rng, functionArgs...)
        end
        parents[i] = population[index]
    end

    return parents
end # function

"""
    selectParentsWithoutReplacement(selector::SelectionOperator,
                                    population::Array{Individual},
                                    evaluator::Evaluator,
                                    rng::Random.AbstractRNG)::Array{Individual}

Performs a selection between individuals without replacement using the comparison
function of the `Evaluator`.
"""
function selectParentsWithoutReplacement(selector::SelectionOperator,
                                         population::Array{Individual},
                                         evaluator::Evaluator,
                                         rng::Random.AbstractRNG)::Array{Individual}

    nSelectedParents = getNSelectedParents(selector)
    compare = getCompareFunction(evaluator)
    popSize = length(population)

    copyPopulation = copyPop(population)
    indexes = collect(1:popSize)

    method = getMethod(selector)
    functionArgs = getFunctionArgs(selector)

    for i=1:nSelectedParents
        if needsComparison(selector)
            index = method(copyPopulation, compare, rng, functionArgs...)
        else
            index = method(copyPopulation, rng, functionArgs...)
        end
        copyPopulation[index] = copyPopulation[end]
        pop!(copyPopulation)
        indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
    end

    parents = population[indexes[end-nSelectedParents+1:end]]

    return parents
end # function

"""
    selectParents_(selector::SelectionOperator, population::Array{Individual},
                   evaluator::Evaluator, rng::Random.AbstractRNG)

Performs the selection of individuals given by the selection method and using
the comparison function between individuals of the `Evaluator`.
"""
function selectParents_(selector::SelectionOperator, population::Array{Individual},
                        evaluator::Evaluator, rng::Random.AbstractRNG)

    if samplingWithReplacement(selector)
        selectParentsWithReplacement(selector, population, evaluator, rng)
    else
        selectParentsWithoutReplacement(selector, population, evaluator, rng)
    end
end # function

"""
    randomSelector(population::Array{Individual}, rng::Random.AbstractRNG)

Performs a random selection between the individuals of a population.

# Mandatory Arguments
- `population::Array{Individual}`: population of individuals that are going to
    be selected.
- `rng::Random.AbstractRNG`: random number generator of the experiment.

# Returns
The index of the selected individual.

See also: [`setSelector`](@ref)
"""
function randomSelector(population::Array{Individual}, rng::Random.AbstractRNG)
    rand(rng, UInt16)%(length(population)) + 1
end # function

"""
    tournamentSelector(population::Array{Individual}, compare::Function,
                       rng::Random.AbstractRNG, k::Integer)

Performs a tournament selection between the individuals of a population.

# Mandatory Arguments
- `population::Array{Individual}`: population of individuals that are going to
    be selected.
- `compare::Function`: `Evaluator`'s comparison method between individuals to
    decide which ones win the tournament.
- `rng::Random.AbstractRNG`: random number generator of the experiment.

# Arguments
- `k::Integer`: number of individuals that compete in the tournament. If `k`
    equals 1, it will be a random selection between the individuals (see
    [`randomSelector`](@ref)).

# Returns
The index of the selected individual.

See also: [`setSelector`](@ref)
"""
function tournamentSelector(population::Array{Individual}, compare::Function,
                            rng::Random.AbstractRNG, k::Integer)
    nInd = length(population) # Number of individuals of the population
    bestIndex = rand(rng, UInt16) % nInd + 1
    best = population[bestIndex]

    for i=2:k
        currentIndex = rand(rng, UInt16) % nInd + 1
        current = population[currentIndex]
        bestIndex = compare(best, current) ? bestIndex : currentIndex
        best = population[bestIndex]
    end

    return bestIndex
end # function

"""
function tournamentSelector(fitnesses::Array, compare::Function, rng::Random.AbstractRNG, k::Integer)
    nInd = length(fitnesses) # Number of individuals of the population
    bestIndex = rand(rng, UInt) % nInd + 1
    best = fitnesses[bestIndex]

    for i=2:k
        currentIndex = rand(rng, UInt) % nInd + 1
        current = fitnesses[currentIndex]
        bestIndex = compare(best, current) ? bestIndex : currentIndex
        best = fitnesses[bestIndex]
    end

    return bestIndex
end
"""

function rouletteSelector(fitnesses::Array, compare::Function, rng::Random.AbstractRNG, k::Integer=-1)

    minimize = compare(1,2)
    popSize = length(fitnesses)
    if minimize
        max = maximum(fitnesses)
        fitnesses = map(x->abs(x-max), fitnesses)
    end

    if k < 0
        k = popSize
        indexes = collect(1:popSize)
    else
        if k < 50
            indexes = randomIndexSelection(popSize, k, rng)
        else
            indexes = randomIndexSelection2(popSize, k, rng)
        end
        fitnesses = fitnesses[indexes]
    end

    sumFitness = sum(fitnesses)

    selected = 1
    random = rand(rng, 0:sumFitness)
    random -= fitnesses[1]

    while random >= 0
        selected += 1
        random -= fitnesses[selected]
    end

    return indexes[selected]
end

function rankingSelector(fitnesses::Array, compare::Function, rng::Random.AbstractRNG, k::Integer=-1)

    minimize = compare(1,2)
    popSize = length(fitnesses)
    if minimize
        newFitnesses = copy(fitnesses)
        max = maximum(newFitnesses)
        newFitnesses = map(x->abs(x-max), newFitnesses)
    end

    if k < 0
        k = popSize
        indexes = collect(1:popSize)
        newFitnesses = fitnesses
    else
        if k < 50
            indexes = randomIndexSelection(popSize, k, rng)
        else
            indexes = randomIndexSelection2(popSize, k, rng)
        end
        newFitnesses = fitnesses[indexes]
    end

    sumFitness = sum(fitnesses)

    selected = 1
    random = rand(rng, 0:sumFitness)
    random -= fitnesses[1]

    while random >= 0
        selected += 1
        random -= fitnesses[selected]
    end

    return indexes[selected]
end

"""
function selectParentsWithoutReplacement(selector::SelectionOperator,
                                         population::Array{Individual},
                                         nParents::Integer, evaluator::Evaluator, rng::Random.AbstractRNG)

    nSelectedParents = getNSelectedParents(selector)
    compare = getCompareFitness(evaluator)
    popSize = length(population)

    if !(typeof(nSelectedParents) <: Integer)
        nSelectedParents = convert(Integer, round(nSelectedParents * popSize))
        remainder = nSelectedParents % nParents

        if remainder != 0
            nSelectedParents = nSelectedParents + nParents - remainder
        end

        if nSelectedParents > popSize
            nSelectedParents -= nParents
        end
    end

    indexes = collect(1:popSize)

    if selector._individualMode == "fitness"
        fitnesses = getFitness(population)

        for i=1:nSelectedParents
            if needsComparison(selector)
                index = selector._method(fitnesses, compare, rng, selector._varArgs...)
                fitnesses[index] = fitnesses[end]
                pop!(fitnesses)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            else
                index = selector._method(fitnesses, rng, selector._varArgs...)
                fitnesses[index] = fitnesses[end]
                pop!(fitnesses)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            end
        end
    elseif selector._individualMode == "representation"
        indReps = getGenotype(population)
        for i=1:nSelectedParents
            if needsComparison(selector)
                index = selector._method(indReps, compare, rng, selector._varArgs...)
                indReps[index] = indReps[end]
                pop!(indReps)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            else
                index = selector._method(indReps, rng, selector._varArgs...)
                indReps[index] = indReps[end]
                pop!(indReps)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            end
        end
    else
        fitnesses = getFitness(population)
        indReps = getGenotype(population)
        for i=1:nSelectedParents
            if needsComparison(selector)
                index = selector._method(indReps, fitnesses, compare, rng, selector._varArgs...)
                indReps[index] = indReps[end]
                pop!(indReps)
                fitnesses[index] = fitnesses[end]
                pop!(fitnesses)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            else
                index = selector._method(indReps, fitnesses, rng, selector._varArgs...)
                indReps[index] = indReps[end]
                pop!(indReps)
                fitnesses[index] = fitnesses[end]
                pop!(fitnesses)
                indexes[index], indexes[end-i+1] = indexes[end-i+1], indexes[index]
            end
        end
    end

    parents = deepcopy(population[indexes[end-nParents+1:end]])

    return parents
end # function
"""

max_stars_repo_name: SergioGmezM/GeneticJulia

max_issues_repo_name: SergioGmezM/GeneticJulia

max_forks_repo_name: SergioGmezM/GeneticJulia


text:
using MITgcmTools, MeshArrays, Plots

p="./hs94.cs-32x32x5/run/"
pp=tempdir()*"/"

(γ,Γ)=GridOfOnes("CubeSphere",6,32)
#f2d(x)=γ.read(permutedims(Float64.(x),[2,1]),MeshArray(γ))
#f3d(x)=γ.read(permutedims(Float64.(x),[2,1,3]),MeshArray(γ))

list_n=("XC","XG","YC","YG","RAC","RAW","RAS","RAZ","DXC","DXG","DYC","DYG","Depth");
for n in list_n
#    tmp=f2d(read_mdsio(p,n))
#    [Γ[n][i].=tmp[i] for i in eachindex(tmp)]
    tmp=read_mdsio(p,n)
    [Γ[n][i].=tmp[(1:32).+(i-1)*32,:] for i in 1:6]
end

#XC=f2d(read_mdsio(p,"XC"))
#YC=f2d(read_mdsio(p,"YC"))
#[Γ["XC"][i].=XC[i] for i in eachindex(XC)];
#[Γ["YC"][i].=YC[i] for i in eachindex(YC)];

lon=[i for i=-179.5:1.0:179.5, j=-89.5:1.0:89.5]
lat=[j for i=-179.5:1.0:179.5, j=-89.5:1.0:89.5]
(f,i,j,w,j_f,j_x,j_y)=InterpolationFactors(Γ,vec(lon),vec(lat))

##

ff=readdir(p); fil="T.0000"
ff=ff[findall(occursin.(fil,ff).*occursin.(".data",ff))]

T=MeshArray(γ)

function myplot(fil)
    tmp=read_mdsio(p,fil)[:,:,3]
    [T[i].=tmp[(1:32).+(i-1)*32,:] for i in 1:6]
    TT=Interpolate(T,f,i,j,w)
    contourf(vec(lon[:,1]),vec(lat[1,:]),TT,clims=(260.,320.))
end

#myplot("T.0000000000")

nt=length(ff)
#nt=300
#dt=Int(ceil(nt/100))
dt=1
anim = @animate for i ∈ 1:dt:nt
    myplot(ff[i])
end

gif(anim,pp*"hs94.cs.10days.gif", fps = 8)

ext: jl


text:
"""
Given substrate and product stoichiometric matrices, rate-expressions and lists
of species and parameters, return a ReactionSystem that describes the chemical
reaction network.

Assumed that the substrate and product stoichiometry matrices are stored as
numspecies by numrxs matrices, with entry (i,j) giving the stoichiometric
coefficient of species i within rx j.
"""
struct MatrixNetwork{S ,T ,U ,V ,W,X} <: NetworkFileFormat
    """The symbolic expressions for each reaction rate."""
    rateexprs::S

    """substoich[i,j] = is the substrate stoichiometric coefficient matrix"""
    substoich::T

    """prodstoich[i,j] is the product stoichiometric coefficient matrix"""
    prodstoich::U

    """species in the network """
    species::V

    """ Parameters """
    params::W

    """independent variable, time """
    t::X
end
MatrixNetwork(rateexprs,substoich,prodstoich; species=Any[],params=Any[],t=nothing) =
           MatrixNetwork(rateexprs,substoich,prodstoich,species,params,t)

# for dense matrices
function loadrxnetwork(mn::MatrixNetwork{S,T,U,V,W,X}; name = gensym(:ReactionSystem)) where {S <: AbstractVector,
        T <: Matrix, U <: Matrix{Int}, V <: AbstractVector, W <: AbstractVector,X <: Any}

    sz = size(mn.substoich)
    @assert sz == size(mn.prodstoich)
    numspecs= sz[1]
    numrxs = sz[2]

    t = (mn.t === nothing) ? Catalyst.DEFAULT_IV : mn.t
    species = isempty(mn.species) ? [funcsym(:S,t,i) for i = 1:numspecs] : mn.species

    # create the reactions
    # we need to create new vectors each time as the ReactionSystem
    # takes ownership of them
    rxs = Vector{Reaction}(undef,numrxs)
    for j = 1:numrxs
        subs    = Any[]
        sstoich = Vector{eltype(mn.substoich)}()
        prods   = Any[]
        pstoich = Vector{eltype(mn.prodstoich)}()

        # stoich
        for i = 1:numspecs
            scoef = mn.substoich[i,j]
            if (scoef > zero(scoef))
                push!(subs, species[i])
                push!(sstoich, scoef)
            end

            pcoef = mn.prodstoich[i,j]
            if (pcoef > zero(pcoef))
                push!(prods, species[i])
                push!(pstoich, pcoef)
            end
        end

        rxs[j] = Reaction(mn.rateexprs[j], subs, prods, sstoich, pstoich)
    end

    ParsedReactionNetwork(ReactionSystem(rxs, t, species, mn.params; name=name))
end

# for sparse matrices
function loadrxnetwork(mn::MatrixNetwork{S,T,U,V,W,X}; name = gensym(:ReactionSystem)) where {S<:AbstractVector,
        T<:SparseMatrixCSC,U<:SparseMatrixCSC{Int,Int},V<:AbstractVector, W<:AbstractVector,X <: Any}
    sz = size(mn.substoich)
    @assert sz == size(mn.prodstoich)
    numspecs= sz[1]
    numrxs = sz[2]

    t       = (mn.t === nothing) ? Catalyst.DEFAULT_IV : mn.t
    species = isempty(mn.species) ? [funcsym(:S,t,i) for i = 1:numspecs] : mn.species

    # create the reactions
    srows = rowvals(mn.substoich)
    svals = nonzeros(mn.substoich)
    prows = rowvals(mn.prodstoich)
    pvals = nonzeros(mn.prodstoich)
    rxs   = Vector{Reaction}(undef,numrxs)
    for j = 1:numrxs
        subs    = Any[]
        sstoich = Vector{eltype(mn.substoich)}()
        prods   = Any[]
        pstoich = Vector{eltype(mn.prodstoich)}()

        for ir in nzrange(mn.substoich, j)
           i     = srows[ir]
           scoef = svals[ir]
           if scoef > zero(scoef)
                push!(subs, species[i])
                push!(sstoich, scoef)
           end
        end

        for ir in nzrange(mn.prodstoich, j)
            i     = prows[ir]
            pcoef = pvals[ir]
            if pcoef > zero(pcoef)
                push!(prods, species[i])
                push!(pstoich, pcoef)
            end
        end

        rxs[j] = Reaction(mn.rateexprs[j], subs, prods, sstoich, pstoich)
     end

     ParsedReactionNetwork(ReactionSystem(rxs, t, species, mn.params; name=name))
end

"""
Given complex stoichiometric matrix ,incidence matrix ,rate-expressions and list of species,parameters
, return a ReactionSystem that describes the chemical reaction network

Notes:
- The column of complex stoichiometric represents composition of reaction complexes,
  with positive entries of size num_of_species by num_of_complexes, where
  the non-zero positive entries in the k'th column denote stoichiometric
  coefficients of the species participating in the k'th reaction complex.

- The complex incidence matrix, is number of complexes by number of reactions with
  Bᵢⱼ = -1, if the i'th complex is the substrate of the j'th reaction,
         1, if the i'th complex is the product of the j'th reaction,
         0, otherwise
"""
struct ComplexMatrixNetwork{S ,T ,U ,V ,W,X} <: NetworkFileFormat
    """The symbolic expressions for each reaction rate."""
    rateexprs::S

    """stoichmat[i,j] = is the stoichiometric coefficient in the j'th reaction for the i'th species"""
    stoichmat::T

    """incidencemat[i,j] is the incidence matrix with incidencemat[i,j] = -1, if the i'th complex
    is the substrate of the j'th reaction, 1, if the i'th complex is the product
    of the j'th reaction, 0, otherwise"""
    incidencemat::U  # all elements are integer always in incidence matrix

    """species in the network """
    species::V

    """ Parameters """
    params::W

    """independent variable, time """
    t::X
end
ComplexMatrixNetwork(rateexprs,stoichmat,incidencemat; species=Any[],params=Any[],t=nothing) =
           ComplexMatrixNetwork(rateexprs,stoichmat,incidencemat,species,params,t)

# for Dense matrices version
function loadrxnetwork(cmn::ComplexMatrixNetwork{S,T,U,V,W,X}; name = gensym(:ReactionSystem)) where {S <: AbstractVector,
                T <: Matrix, U <: Matrix{Int}, V <: AbstractVector,
                                        W <: AbstractVector,X <: Any}
    numspecs, numcomp = size(cmn.stoichmat)
    @assert all(>=(0),cmn.stoichmat)
    @assert numcomp == size(cmn.incidencemat, 1)
    @assert all(∈([-1,0,1]),cmn.incidencemat)
    numrxs = size(cmn.incidencemat, 2)

    t = (cmn.t === nothing) ? Catalyst.DEFAULT_IV : cmn.t
    species = isempty(cmn.species) ? [funcsym(:S,t,i) for i = 1:numspecs] : cmn.species

    rxs = Vector{Reaction}(undef,numrxs)
    sc_ind = argmin(cmn.incidencemat, dims=1)  # cartesian indices of substrate complexes
    pc_ind = argmax(cmn.incidencemat, dims=1)  # cartesian indicies of products complexes
    for i ∈ 1:numrxs

        # substrate index for i'th reaction in species(rn)
        ss_ind = findall(!iszero, @view cmn.stoichmat[:,sc_ind[i][1]])

        # products index for i'th reaction in species(rn)
        ps_ind = findall(!iszero, @view cmn.stoichmat[:,pc_ind[i][1]])

        if isempty(ss_ind) && !isempty(ps_ind)
            rxs[i] = Reaction(cmn.rateexprs[i], nothing, species[ps_ind],
                        nothing, cmn.stoichmat[ps_ind,pc_ind[i][1]])

        elseif !isempty(ss_ind) &&  isempty(ps_ind)
            rxs[i] = Reaction(cmn.rateexprs[i], species[ss_ind], nothing,
                        cmn.stoichmat[ss_ind,sc_ind[i][1]], nothing)
        else
            rxs[i] = Reaction(cmn.rateexprs[i], species[ss_ind],species[ps_ind],
                        cmn.stoichmat[ss_ind,sc_ind[i][1]],
                        cmn.stoichmat[ps_ind,pc_ind[i][1]])
        end
    end

    ParsedReactionNetwork(ReactionSystem(rxs, t, species, cmn.params; name=name))
end

# for sparse matrices version
function loadrxnetwork(cmn::ComplexMatrixNetwork{S,T,U,V,W,X}; name = gensym(:ReactionSystem)) where {S<:AbstractVector,
                T<:SparseMatrixCSC,U<:SparseMatrixCSC{Int,Int},V<:AbstractVector,
                                                    W<:AbstractVector,X <: Any}
    numspecs, numcomp = size(cmn.stoichmat)
    @assert all(>=(0),cmn.stoichmat)
    @assert numcomp == size(cmn.incidencemat, 1)
    @assert all(∈([-1,0,1]),cmn.incidencemat)
    numrxs = size(cmn.incidencemat, 2)

    t = (cmn.t === nothing) ? Catalyst.DEFAULT_IV : cmn.t
    species = isempty(cmn.species) ? [funcsym(:S,t,i) for i = 1:numspecs] : cmn.species

    rxs = Vector{Reaction}(undef,numrxs)
    sc_ind = argmin(cmn.incidencemat, dims=1)  # cartesian indices of substrate complexes
    pc_ind = argmax(cmn.incidencemat, dims=1)  # cartesian indicies of products complexes
    rows = rowvals(cmn.stoichmat)
    vals = nonzeros(cmn.stoichmat)
    for i ∈ 1:numrxs
        # substrate index for i'th reaction in species(rn)
        ss_ind = @view rows[nzrange(cmn.stoichmat, sc_ind[i][1])]
        # products index for i'th reaction in species(rn)
        ps_ind = @view rows[nzrange(cmn.stoichmat, pc_ind[i][1])]

        if isempty(ss_ind) && !isempty(ps_ind)
            rxs[i] = Reaction(cmn.rateexprs[i], nothing, species[ps_ind],
                        nothing, vals[nzrange(cmn.stoichmat, pc_ind[i][1])])

        elseif !isempty(ss_ind) &&  isempty(ps_ind)
            rxs[i] = Reaction(cmn.rateexprs[i], species[ss_ind], nothing,
                        vals[nzrange(cmn.stoichmat, sc_ind[i][1])], nothing)
        else
            rxs[i] = Reaction(cmn.rateexprs[i], species[ss_ind],species[ps_ind],
                        vals[nzrange(cmn.stoichmat, sc_ind[i][1])],
                        vals[nzrange(cmn.stoichmat, pc_ind[i][1])])
        end
    end

    ParsedReactionNetwork(ReactionSystem(rxs, t, species, cmn.params; name=name))
end


text:
@inline function getindex(sa::HybridArray{S}, ::Colon) where S
    return HybridArray{S}(getindex(sa.data, :))
end

Base.@propagate_inbounds function getindex(sa::HybridArray{S}, inds::Int...) where S
    return getindex(sa.data, inds...)
end

Base.@propagate_inbounds function getindex(sa::HybridArray{S}, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where S
    _getindex(all_dynamic_fixed_val(S, inds...), sa, inds...)
end

Base.@propagate_inbounds function _getindex(::Val{:dynamic_fixed_true}, sa::HybridArray, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...)
    return _getindex_all_static(sa, inds...)
end

function _get_indices(i::Tuple{}, j::Int)
    return ()
end

function _get_indices(i::Tuple, j::Int, i1::Type{Int}, inds...)
    return (:(inds[$j]), _get_indices(i, j+1, inds...)...)
end

function _get_indices(i::Tuple, j::Int, i1::Type{T}, inds...) where T<:StaticArray{<:Tuple, Int}
    return (:(inds[$j][$(i[1])]), _get_indices(i[2:end], j+1, inds...)...)
end

function _get_indices(i::Tuple, j::Int, i1::Type{Colon}, inds...)
    return (i[1], _get_indices(i[2:end], j+1, inds...)...)
end

_totally_linear() = true
_totally_linear(inds...) = false
_totally_linear(inds::Type{Int}...) = true
_totally_linear(inds::Type{Colon}...) = true
_totally_linear(i1::Type{Colon}, inds...) = _totally_linear(inds...)

function new_out_size_nongen(::Type{Size}, inds...) where Size
    os = []
    map(Size.parameters, inds) do s, i
        if i == Int
        elseif i <: StaticVector
            push!(os, length(i))
        elseif i == Colon
            push!(os, s)
        else
            error("Unknown index type: $i")
        end
    end
    return tuple(os...)
end

"""
    _get_linear_inds(S, inds...)

Returns linear indices for given size and access indices.
Order is selected to make setindex! with array input be linearly indexed by
position of index in returned vector.
"""
function _get_linear_inds(S, inds...)
    newsize = new_out_size_nongen(S, inds...)
    indices = CartesianIndices(newsize)
    out_inds = Any[]
    sizeprods = Union{Int, Expr}[1]
    needs_dynsize = false
    for (i, s) in enumerate(S.parameters[1:(end-1)])
        if isa(s, Int) && isa(sizeprods[end], Int)
            push!(sizeprods, s*sizeprods[end])
        elseif isa(s, Int)
            push!(sizeprods, :($s*$(sizeprods[end])))
        elseif isa(s, Dynamic)
            needs_dynsize = true
            push!(sizeprods, :(dynsize[$i]*$(sizeprods[end])))
        end
    end

    needs_sizeprod = false
    if _totally_linear(inds...)
        i1 = 0
        for (iik, ik) ∈ enumerate(inds)
            if isa(ik, Type{Int})
                needs_sizeprod = true
                i1 = :($i1 + (inds[$iik]-1)*sizeprods[$iik])
            end
        end
        out_inds = tuple((:(i1 + $k) for k ∈ 1:StaticArrays.tuple_prod(newsize))...)

        preamble = if needs_sizeprod && needs_dynsize
            quote
                dynsize = size(sa)
                sizeprods = tuple($(sizeprods...))
                i1 = $i1
            end
        elseif needs_sizeprod
            quote
                sizeprods = tuple($(sizeprods...))
                i1 = $i1
            end
        else
            quote
                i1 = $i1
            end
        end
        return preamble, out_inds
    else
        return nothing
    end
end

@generated function _getindex_all_static(sa::HybridArray{S,T}, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where {S,T}
    newsize = new_out_size_nongen(S, inds...)
    exprs = Vector{Expr}(undef, length(newsize))

    indices = CartesianIndices(newsize)
    exprs = similar(indices, Expr)
    Tnewsize = Tuple{newsize...}

    lininds = _get_linear_inds(S, inds...)
    if lininds === nothing
        for current_ind ∈ indices
            cinds = _get_indices(current_ind.I, 1, inds...)
            exprs[current_ind.I...] = :(getindex(sadata, $(cinds...)))
        end
        return quote
            Base.@_propagate_inbounds_meta
            sadata = sa.data
            SArray{$Tnewsize,$T}(tuple($(exprs...)))
        end
    else
        exprs = [:(getindex(sadata, $id)) for id ∈ lininds[2]]
        return quote
            Base.@_propagate_inbounds_meta
            sadata = sa.data
            $(lininds[1])
            SArray{$Tnewsize,$T}(tuple($(exprs...)))
        end
    end
end

function new_out_size(S::Type{Size}, inds::StaticArrays.StaticIndexing...) where Size
    return new_out_size(S, map(StaticArrays.unwrap, inds)...)
end

# _get_static_vector_length is used in a generated function so using a generic function
# may not be a good idea
_get_static_vector_length(::Type{<:StaticVector{N}}) where {N} = N

@generated function new_out_size(::Type{Size}, inds...) where Size
    os = []
    map(Size.parameters, inds) do s, i
        if i == Int
        elseif i <: StaticVector
            push!(os, _get_static_vector_length(i))
        elseif i == Colon || i <: Base.Slice
            push!(os, s)
        elseif i <: SOneTo
            push!(os, i.parameters[1])
        elseif i <: Base.OneTo || i <: AbstractArray
            push!(os, Dynamic())
        else
            error("Unknown index type: $i")
        end
    end
    return Tuple{os...}
end

@inline function _getindex(::Val{:dynamic_fixed_false}, sa::HybridArray{S}, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where S
    newsize = new_out_size(S, inds...)
    return HybridArray{newsize}(getindex(sa.data, inds...))
end

# setindex stuff

Base.@propagate_inbounds function setindex!(a::HybridArray, value, inds::Int...)
    Base.@boundscheck checkbounds(a, inds...)
    _setindex!_scalar(Size(a), a, value, inds...)
end

@generated function _setindex!_scalar(::Size{S}, a::HybridArray, value, inds::Int...) where S
    if length(inds) == 0
        return quote
            Base.@_propagate_inbounds_meta
            a[1] = value
        end
    end

    stride = :(1)
    ind_expr = :()
    for i ∈ 1:length(inds)
        if i == 1
            ind_expr = :(inds[1])
        else
            ind_expr = :($ind_expr + $stride * (inds[$i] - 1))
        end
        # it's possible that one of the trailing indices is an additional 1.
        if length(S) < i || isa(S[i], StaticArrays.Dynamic)
            stride = :($stride * size(a.data, $i))
        else
            stride = :($stride * $(S[i]))
        end
    end
    return quote
        Base.@_inline_meta
        Base.@_propagate_inbounds_meta
        a.data[$ind_expr] = value
    end
end

Base.@propagate_inbounds function setindex!(sa::HybridArray{S}, value, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where S
    _setindex!(all_dynamic_fixed_val(S, inds...), sa, value, inds...)
end

@inline function _setindex!(::Val{:dynamic_fixed_false}, sa::HybridArray{S}, value, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where S
    newsize = new_out_size(S, inds...)
    return HybridArray{newsize}(setindex!(sa.data, value, inds...))
end

Base.@propagate_inbounds function _setindex!(::Val{:dynamic_fixed_true}, sa::HybridArray, value, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...)
    return _setindex!_all_static(sa, value, inds...)
end

@generated function _setindex!_all_static(sa::HybridArray{S,T}, v::AbstractArray, inds::Union{Int, StaticArray{<:Tuple, Int}, Colon}...) where {S,T}
    newsize = new_out_size_nongen(S, inds...)
    exprs = Vector{Expr}(undef, length(newsize))

    indices = CartesianIndices(newsize)
    Tnewsize = Tuple{newsize...}

    if v <: StaticArray
        newlen = StaticArrays.tuple_prod(newsize)
        if Length(v) != newlen
            return quote
                throw(DimensionMismatch("tried to assign $(length(v))-element array to $($newlen) destination"))
            end
        end
    end

    lininds = _get_linear_inds(S, inds...)
    if lininds === nothing
        exprs = similar(indices, Expr)
        for current_ind ∈ indices
            cinds = _get_indices(current_ind.I, 1, inds...)
            exprs[current_ind.I...] = :(setindex!(sadata, v[$(current_ind.I...)], $(cinds...)))
        end

        if v <: StaticArray
            return quote
                Base.@_propagate_inbounds_meta
                sadata = sa.data
                @inbounds $(Expr(:block, exprs...))
            end
        else
            return quote
                Base.@_propagate_inbounds_meta
                if size(v) != $newsize
                    throw(DimensionMismatch("tried to assign array of size $(size(v)) to destination of size $($newsize)"))
                end
                sadata = sa.data
                @inbounds $(Expr(:block, exprs...))
            end
        end

    else
        exprs = [:(setindex!(sadata, v[$iid], $id)) for (iid, id) ∈ enumerate(lininds[2])]

        if v <: StaticArray
            return quote
                Base.@_propagate_inbounds_meta
                $(lininds[1])
                sadata = sa.data
                @inbounds $(Expr(:block, exprs...))
            end
        else
            return quote
                Base.@_propagate_inbounds_meta
                if size(v) != $newsize
                    throw(DimensionMismatch("tried to assign array of size $(size(v)) to destination of size $($newsize)"))
                end
                $(lininds[1])
                sadata = sa.data
                @inbounds $(Expr(:block, exprs...))
            end
        end
    end
end

@inline function _view_hybrid(a::HybridArray{S}, ::Val{:dynamic_fixed_true}, inner_view, indices...) where {S}
    new_size = new_out_size(S, indices...)
    return SizedArray{new_size}(inner_view)
end

@inline function _view_hybrid(a::HybridArray{S}, ::Val{:dynamic_fixed_false}, inner_view, indices...) where {S}
    new_size = new_out_size(S, indices...)
    return HybridArray{new_size}(inner_view)
end

@inline function Base.view(
    a::HybridArray{S},
    indices...,
) where {S}
    inner_view = invoke(view, Tuple{AbstractArray, typeof(indices).parameters...}, a, indices...)
    return _view_hybrid(a, all_dynamic_fixed_val(S, indices...), inner_view, indices...)
end


text:
using COSMO, Test, Random

settings = COSMO.Settings()
iter = 10
cost = 20.
r_prim = 1.5e-3
r_dual = 1.2e-2
status = :Solved
rt = 0.7

@testset "Printing" begin
   @test COSMO.print_iteration(settings, iter, cost, r_prim, r_dual) == nothing
   @test COSMO.print_iteration(settings, settings.check_termination, cost, r_prim, r_dual) == nothing
   @test COSMO.print_result(status, iter, cost, rt) == nothing
end
nothing


text:
# Set the listener properly
using ZMQ
using MAT

ctx=Context()
s1=Socket(ctx, REP)

ZMQ.bind(s1, "ipc:///tmp/armopt")
println("Waiting for a message")
msg = ZMQ.recv(s1)
#out=convert(IOStream, msg)
#seek(out,0)
#read out::MemIO as usual, eg. read(out,...) or takebuf_string(out)
#or, conveniently, use bytestring(msg) to retrieve a string

println("Received a message!")

matfile = bytestring(msg)
println("I received ", matfile)

file = matopen(matfile)
varnames = names(file)
println(varnames)
close(file)

# This works fine
vars = matread(matfile)
#println(vars)

# Original Path
qw0 = vars["qwPath"]
# Task Space target
vw0 = vars["vwPath"]
# Task Space effective
vw0 = vars["vwEffective"]

ZMQ.send(s1, Message(matfile))
ZMQ.close(s1)

ZMQ.close(ctx)

#
# using ECOS, SCS, Gurobi
# using Convex
#
# mysolver = ECOSSolver()
# #mysolver = SCSSolver()
# #mysolver = GurobiSolver()
#
# items = [:Gold, :Silver, :Bronze]
# values = [5.0, 3.0, 1.0]
# weights = [2.0, 1.5, 0.3]
#
# # Define a variable of size 3, each index representing an item
# x = Convex.Variable(3)
# p = maximize(x' * values, 0 <= x, x <= 1, x' * weights <= 3)
# @time Convex.solve!(p, mysolver)
#
# println("Optimal objective: ", p.optval)
# items = Dict{Symbol,Float64}(zip(items, x.value))
#
# # [:Gold 0.9999971880377178
# #  :Silver 0.46667637765641057
# #  :Bronze 0.9999998036351865]

max_stars_repo_name: ToyotaResearchInstitute/rad-robot

max_issues_repo_name: ToyotaResearchInstitute/rad-robot

max_forks_repo_name: ToyotaResearchInstitute/rad-robot


text:
using BinaryBuilder

name = "Notcurses"
version = v"3.0.2"
sources = [
    GitSource("https://github.com/dankamongmen/notcurses",
              "c19a2477ff41586bbe899f354d7d65b54445be2a"),
    DirectorySource("./bundled"),
]

script = raw"""
cd ${WORKSPACE}/srcdir/notcurses*/
atomic_patch -p1 ${WORKSPACE}/srcdir/patches/repent.patch
atomic_patch -p1 ${WORKSPACE}/srcdir/patches/0001-also-look-for-shared-libraries-on-Windows.patch

if [[ $target == *mingw* ]]; then
    export CFLAGS="${CFLAGS} -D_WIN32_WINNT=0x0600"
    cp ${WORKSPACE}/srcdir/headers/pthread_time.h "/opt/${target}/${target}/sys-root/include/pthread_time.h"
fi

install_license COPYRIGHT

mkdir build && cd build

export CFLAGS="${CFLAGS} -I${includedir}"

FLAGS=(-DCMAKE_TOOLCHAIN_FILE=${CMAKE_TARGET_TOOLCHAIN}
       -DCMAKE_INSTALL_PREFIX=${prefix}
       -DCMAKE_BUILD_TYPE=Release
       -DBUILD_SHARED_LIBS=ON
       -DUSE_DOCTEST=off
       -DUSE_PANDOC=off
       -DUSE_STATIC=off
       -DUSE_QRCODEGEN=off
       -DBUILD_EXECUTABLES=off
       -DUSE_POC=off
       -DUSE_MULTIMEDIA=none
       )

cmake .. "${FLAGS[@]}"
make -j${nproc}
make install
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line.
platforms = supported_platforms(; experimental=true)

# The products that we will ensure are always built.
products = [
    LibraryProduct("libnotcurses", :libnotcurses)
    LibraryProduct("libnotcurses-core", :libnotcurses_core)
    LibraryProduct("libnotcurses-ffi", :libnotcurses_ffi)
]

# Dependencies that must be installed before this package can be built.
dependencies = [
    Dependency("Ncurses_jll"),
    Dependency("libunistring_jll"),
    Dependency("libdeflate_jll"),
]

build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; preferred_gcc_version=v"7", julia_compat="1.6")

max_stars_repo_name: TeroFrondelius/Yggdrasil

max_issues_repo_name: TeroFrondelius/Yggdrasil

max_forks_repo_name: TeroFrondelius/Yggdrasil


text:
using Documenter, Augmentor, Random

# Autogenerate documentation markdown and jupyter notebooks
# for all the scripts in the "examples/" subfolder.
include("exampleweaver.jl")
#ExampleWeaver.weave(overwrite=false, execute=true)
ExampleWeaver.weave(overwrite=true, execute=true)

# Define the documentation order of the operations. The whole
# purpose of this vector is literally just to dictate in what
# chronological order the operations are documented.
op_fnames = [
    "flipx",
    "flipy",
    "rotate90",
    "rotate270",
    "rotate180",
    "rotate",
    "shearx",
    "sheary",
    "scale",
    "zoom",
    "elasticdistortion",
    "crop",
    "cropnative",
    "cropsize",
    "cropratio",
    "rcropratio",
    "resize",
    "converteltype",
    "mapfun",
    "aggmapfun",
    "splitchannels",
    "combinechannels",
    "permutedims",
    "reshape",
    "noop",
    "cacheimage",
    "either",
]
dict_order = Dict(fname * ".md" => i for (i, fname) in enumerate(op_fnames))
myless(a, b) = dict_order[a] < dict_order[b]

# --------------------------------------------------------------------

Random.seed!(1337)
makedocs(
    modules = [Augmentor],
    clean = false,
    format = :html,
    assets = [
        joinpath("assets", "favicon.ico"),
        joinpath("assets", "style.css"),
    ],
    sitename = "Augmentor.jl",
    authors = "Christof Stocker",
    linkcheck = !("skiplinks" in ARGS),
    pages = Any[
        "Home" => "index.md",
        "gettingstarted.md",
        "Introduction and Motivation" => Any[
            "background.md",
            "images.md",
        ],
        "User's Guide" => Any[
            "interface.md",
            hide("operations.md", Any[joinpath("operations", fname) for fname in sort(readdir(joinpath(@__DIR__, "src", "operations")), lt = myless) if splitext(fname)[2] == ".md"]),
        ],
        "Tutorials" => joinpath.("generated", ExampleWeaver.listmarkdown()),
        hide("Indices" => "indices.md"),
        "LICENSE.md",
    ],
    html_prettyurls = !("local" in ARGS),
)

#= Remove this in final version
deploydocs(
    repo = "github.com/Evizero/Augmentor.jl.git",
    target = "build",
    julia = "1.0",
    deps = nothing,
    make = nothing,
)

# --------------------------------------------------------------------
# Post-process the generated HTML files of the examples/tutorials
# 1. Redirect "Edit on Github" link to the "examples/*.jl" file
# 2. Add a link in the top right corner to the Juypter notebook

build_dir = abspath(joinpath(@__DIR__, "build"))
for markdownname in ExampleWeaver.listmarkdown()
    name = splitext(markdownname)[1]
    htmlpath = joinpath(build_dir, "generated", name, "index.html")
    str_html = read(htmlpath, String)
    # replace github url to .jl file
    str_html = replace(
        str_html,
        r"docs/src/generated/([^.]*)\.md" =>
        s"examples/\1.jl"
    )
    # insert link to jupyter notebook
    str_html = replace(
        str_html,
        r"(<a class=\"edit-page\".*GitHub<\/a>)" =>
        s"\1  Juypter Notebook"
    )
    href = "https://nbviewer.jupyter.org/github/Evizero/Augmentor.jl/blob/gh-pages/generated/$(name * ".ipynb")"
    str_html = replace(str_html, "___HREFPLACEHOLDER___" => href)
    # overwrite html file
    write(htmlpath, str_html)
end

=#


text:
# Modelo de condensador com duas alimentacoes
type Condensador_2alim
	Condensador_2alim()=begin
		PP2=outers.PP2
		new(
			DanaPlugin(Dict{Symbol,Any}(
				:Brief=>"Steam tables"
			)),
			Potencia(Dict{Symbol,Any}(
				:Brief=>"Taxa de calor removido"
			)),
			Dif_Temp(Dict{Symbol,Any}(
				:Brief=>"Grau de sub-resfriamento"
			)),
			Corrente(Dict{Symbol,Any}(
				:Brief=>"Corrente com pressao igual a saida",
				:Symbol=>"_{in1}"
			)),
			Corrente(Dict{Symbol,Any}(
				:Symbol=>"_{in2}"
			)),
			Corrente(Dict{Symbol,Any}(
				:Symbol=>"_{out}"
			)),
			[
				:(Fout.P = Fin1.P),
				:(Fout.T = PP2.Tsat(Fout.P) - G_S),
				:([Fout.S,Fout.H] = PP2.propPTl(Fout.P,Fout.T)),
				:(Fout.F = Fin1.F + Fin2.F),
				:(Q_COND = Fin1.F * Fin1.H + Fin2.F * Fin2.H - Fout.F * Fout.H),
			],
			[
				"","","","","",
			],
			[:PP2,],
			[:Q_COND,:G_S,:Fin1,:Fin2,:Fout,]
		)
	end
	PP2::DanaPlugin
	Q_COND::Potencia
	G_S::Dif_Temp
	Fin1::Corrente
	Fin2::Corrente
	Fout::Corrente
	equations::Array{Expr,1}
	equationNames::Array{String,1}
	parameters::Array{Symbol,1}
	variables::Array{Symbol,1}
end
export Condensador_2alim
function setEquationFlow(in::Condensador_2alim)
	addEquation(1)
	addEquation(2)
	addEquation(3)
	addEquation(4)
	addEquation(5)
end


text:
@with_kw mutable struct LShapedData{T <: AbstractFloat}
    Q::T = 1e10
    θ::T = -1e10
    master_objective::AffineDecisionFunction{T} = zero(AffineDecisionFunction{T})
    no_objective::Bool = false
    num_cuts::Int = 0
    iterations::Int = 1
    consolidations::Int = 0
end

@with_kw mutable struct LShapedParameters{T <: AbstractFloat}
    τ::T = 1e-6
    cut_scaling::T = 1.0
    debug::Bool = false
    time_limit::T = Inf
    log::Bool = true
    keep::Bool = true
    offset::Int = 0
    indent::Int = 0
end

"""
    LShapedAlgorithm

Functor object for the L-shaped algorithm.

...
# Algorithm parameters
- `τ::AbstractFloat = 1e-6`: Relative tolerance for convergence checks.
- `debug::Bool = false`: Specifies if extra information should be saved for debugging purposes. Defaults to false for memory efficiency.
- `cut_scaling::AbstractFloat = 1.0`: Rescaling factor for cutting planes to improve numerical stability.
- `log::Bool = true`: Specifices if L-shaped procedure should be logged on standard output or not.
...
"""
struct LShapedAlgorithm{T <: AbstractFloat,
                        A <: AbstractVector,
                        ST <: VerticalStructure,
                        M <: MOI.AbstractOptimizer,
                        E <: AbstractLShapedExecution,
                        F <: AbstractFeasibility,
                        R <: AbstractRegularization,
                        Agg <: AbstractAggregation,
                        C <: AbstractConsolidation} <: AbstractLShaped
    structure::ST
    num_subproblems::Int
    data::LShapedData{T}
    parameters::LShapedParameters{T}

    # Master
    master::M
    decisions::Decisions
    x::A
    Q_history::A

    # Execution
    execution::E

    # Feasibility
    feasibility::F

    # Regularization
    regularization::R

    # Cuts
    master_variables::Vector{MOI.VariableIndex}
    cut_constraints::Vector{CutConstraint}
    cuts::Vector{AnySparseOptimalityCut{T}}
    aggregation::Agg
    consolidation::C
    θ_history::A

    progress::ProgressThresh{T}

    function LShapedAlgorithm(structure::VerticalStructure,
                              x₀::AbstractVector,
                              feasibility_cuts::Bool,
                              _execution::AbstractExecution,
                              regularizer::AbstractRegularizer,
                              aggregator::AbstractAggregator,
                              consolidator::AbstractConsolidator; kw...)
        # Sanity checks
        length(x₀) != num_decisions(structure, 1) && error("Incorrect length of starting guess, has ", length(x₀), " should be ", num_decisions(structure, 1))
        n = num_subproblems(structure, 2)
        n == 0 && error("No subproblems in stochastic program. Cannot run L-shaped procedure.")
        # Float types
        T = promote_type(eltype(x₀), Float32)
        x₀_ = convert(AbstractVector{T}, copy(x₀))
        A = typeof(x₀_)
        # Structure
        ST = typeof(structure)
        M = typeof(backend(structure.first_stage))
        # Feasibility
        feasibility = feasibility_cuts ? HandleFeasibility(T) : IgnoreFeasibility()
        F = typeof(feasibility)
        # Execution policy
        execution = _execution(structure, F, T, A)
        E = typeof(execution)
        # Regularization policy
        regularization = regularizer(structure.decisions[1], x₀_)
        R = typeof(regularization)
        # Aggregation policy
        aggregation = aggregator(n, T)
        Agg = typeof(aggregation)
        # Consolidation policy
        consolidation = consolidator(T)
        C = typeof(consolidation)
        # Algorithm parameters
        params = LShapedParameters{T}(; kw...)

        lshaped = new{T,A,ST,M,E,F,R,Agg,C}(structure,
                                            n,
                                            LShapedData{T}(),
                                            params,
                                            backend(structure.first_stage),
                                            structure.decisions[1],
                                            x₀_,
                                            A(),
                                            execution,
                                            feasibility,
                                            regularization,
                                            Vector{MOI.VariableIndex}(),
                                            Vector{CutConstraint}(),
                                            Vector{SparseOptimalityCut{T}}(),
                                            aggregation,
                                            consolidation,
                                            A(),
                                            ProgressThresh(T(1.0), 0.0, "$(indentstr(params.indent))L-Shaped Gap "))
        # Initialize solver
        initialize!(lshaped)
        return lshaped
    end
end

function show(io::IO, lshaped::LShapedAlgorithm)
    println(io, typeof(lshaped).name.name)
    println(io, "State:")
    show(io, lshaped.data)
    println(io, "Parameters:")
    show(io, lshaped.parameters)
end

function show(io::IO, ::MIME"text/plain", lshaped::LShapedAlgorithm)
    show(io, lshaped)
end

function (lshaped::LShapedAlgorithm)()
    # Reset timer
    lshaped.progress.tfirst = lshaped.progress.tlast = time()
    # Start workers (if any)
    start_workers!(lshaped)
    # Start procedure
    while true
        status = iterate!(lshaped)
        if status !== nothing
            close_workers!(lshaped)
            return status
        end
    end
end


text:
# Interface to Pumi's Geometric Modelling Interface

module gmi

using PumiConfig
using apf_types
using gmi_types

# allow types to be accessed as gmi.foo
import gmi_types: ModelEntity, Model, Gmi_set, Gmi_iter

# iteration
import Base: start, next, done, eltype, iteratorsize, iteratoreltype, eltype
using Base: SizeUnknown, HasEltype

const gmi_lib = joinpath(CONFIG_PATHS["PUMI_LIBDIR"], "libgmi")

#------------------------------------------------------------------------------
# gmi_set stuff

"""
  ccall wrapper.  Users should not call this, call the constructor for
  Gmi_set instead.
"""
function _make_set(n::Integer)

  p = ccall( (:gmi_make_set, gmi_lib), Ptr{Void}, (Cint,), n)
  return p
end

function free_set(s::Ptr{Void})

  ccall( (:gmi_free_set, gmi_lib), Void, (Ptr{Void},), s)
end

function free_set(s::Gmi_set)

  free_set(s.p)
end

function Gmi_set(n::Integer)

  p = _make_set(n)
  v = Vector{ModelEntity}(n)

  obj = Gmi_set(p, v)
  finalizer(obj, free_set)

  return obj
end

#------------------------------------------------------------------------------
# gmi_iter

function _gmi_begin(g::Model, dim::Integer)
  p = ccall( (:gmi_begin, gmi_lib), Ptr{Void}, (Model, Cint), g, dim)
  return _Gmi_iter(p)
end

function _gmi_next(g::Model, it::_Gmi_iter)
  ge = ccall( (:gmi_next, gmi_lib), Ptr{Void}, (Model, _Gmi_iter), g, it)
  return ModelEntity(ge)
end

function gmi_end(g::Model, it::_Gmi_iter)
  ccall( (:gmi_end, gmi_lib), Void, (Model, _Gmi_iter), g, it)
end

function gmi_end(obj::Gmi_iter)
  gmi_end(obj.g, obj.it)
end

function Gmi_iter(g::Model, dim::Integer)
  it = _gmi_begin(g, dim)
  obj = Gmi_iter(it, g)

  return obj
end

function start(obj::Gmi_iter)
  return  _gmi_next(obj.g, obj.it)
end

function next(obj::Gmi_iter, state::ModelEntity)

  # the state always leads the current element by 1, so we can do
  # the done check
  next_el = _gmi_next(obj.g, obj.it)
  return (state, next_el)
end

function done(obj::Gmi_iter, state::ModelEntity)

  isdone = state == NullModelEntity
  if isdone
    gmi_end(obj)
  end

  return state == NullModelEntity
end

function iteratorsize(::Type{Gmi_iter})
  return SizeUnknown()
end

function iteratoreltype(::Type{Gmi_iter})
  return HasEltype()
end

function eltype(::Type{Gmi_iter})
  return ModelEntity
end

#------------------------------------------------------------------------------
# more methods

function dim(g::Model, ge::ModelEntity)
  return ccall((:gmi_dim, gmi_lib), Cint, (Model, ModelEntity), g, ge)
end

function tag(g::Model, ge::ModelEntity)
  return ccall((:gmi_tag, gmi_lib), Cint, (Model, ModelEntity), g, ge)
end

function find(g::Model, dim::Integer, tag::Integer)
  ge = ccall( (:gmi_find, gmi_lib), Ptr{Void}, (Model, Cint, Cint), g, dim, tag)
  return ModelEntity(ge)
end

function adjacent(g::Model, ge::ModelEntity, dim::Integer)

  # stuff these functions into libPumiInterface
  n = ccall( (:gmi_adjacent_count, pumi_libname), Cint, (Model, ModelEntity, Cint), g, ge, dim)

  v = Vector{ModelEntity}(n)
  gset = ccall( (:gmi_adjacent_get, pumi_libname), Ptr{Void}, (Ptr{ModelEntity},), v)

  obj = Gmi_set(gset, v)
  finalizer(obj, free_set)

  return obj
end

function can_eval(g::Model)

  val = ccall( (:gmi_can_eval, gmi_lib), Cint, (Model,), g)
  return val != 0
end

"""
  Note the prefix `g`, to avoid conflict with Julia's eval
"""
function geval(g::Model, ge::ModelEntity, p::AbstractVector{Cdouble},
               x::AbstractVector{Cdouble})

  @assert length(p) >= 2
  @assert length(x) == 3
  ccall( (:gmi_eval, gmi_lib), Void, (Model, ModelEntity, Ptr{Cdouble}, Ptr{Cdouble}), g, ge, p, x)

  # make sure any unused components are zeroed out
  for i=3:length(p)
    p[i] = 0
  end

end

function reparam(g::Model, ge_from::ModelEntity,
                 from_p::AbstractVector{Cdouble},
                 ge_to::ModelEntity, to_p::AbstractVector{Cdouble})

  ccall( (:gmi_reparam, gmi_lib), Void, (Model, ModelEntity, Ptr{Cdouble}, ModelEntity, Ptr{Cdouble}), g, ge_from,  from_p, ge_to, to_p)
end

function periodic(g::Model, ge::ModelEntity, dim::Integer)

  val = ccall( (:gmi_periodic, gmi_lib), Cint, (Model, ModelEntity, Cint), g, ge, dim)

  return val != 0
end

function range(g::Model, ge::ModelEntity, dim::Integer, r::AbstractVector{Cdouble})

  @assert length(r) == 2

  ccall( (:gmi_range, gmi_lib), Void, (Model, ModelEntity, Cint, Ptr{Cdouble}), g, ge, dim, r)
end

function closest_point(g::Model, ge::ModelEntity,
                           from::AbstractVector{Cdouble},
                           to::AbstractVector{Cdouble},
                           to_p::AbstractVector{Cdouble})

  @assert length(from) == 3
  @assert length(to) == 3
  @assert length(to_p) >= 2

  ccall( (:gmi_closest_point, gmi_lib), Void, (Model, ModelEntity, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}), g, ge, from, to, to_p)

  # make sure any unused components are zeroed out
  for i=3:length(to_p)
    to_p[i] = 0
  end

end

"""
  For GeomSim models, this only works for geometric faces (not edges, even
  if the model is 2D).
"""
function normal(g::Model, ge::ModelEntity, p::AbstractVector{Cdouble},
                n::AbstractVector{Cdouble})

  ccall( (:gmi_normal, gmi_lib), Void, (Model, ModelEntity, Ptr{Cdouble}, Ptr{Cdouble}), g, ge, p, n)
end

"""
  If `ge` is an edge, only `t0` is populated with dX/dp[1].  If `ge` is
  a face`, then `t1` is also populated with dX/dp[2]
"""
function first_derivative(g::Model, ge::ModelEntity, p::AbstractVector{Cdouble},
                      t0::AbstractVector{Cdouble}, t1::AbstractVector{Cdouble})

  ccall( (:gmi_first_derivative, gmi_lib), Void, (Model, ModelEntity, Ptr{Cdouble}, Ptr{Cdouble}, Ptr{Cdouble}), g, ge, p, t0, t1)

end

function is_point_in_region(g::Model, ge::ModelEntity,
                            point::AbstractVector{Cdouble})

  val = ccall( (:gmi_is_point_in_region, gmi_lib), Cint, (Model, ModelEntity, Ptr{Cdouble}), g, ge, point)

  return val != 0
end

function is_in_closure_of(g::Model, ge1::ModelEntity, ge2::ModelEntity)

  val = ccall( (:gmi_is_in_closure_of, gmi_lib), Cint, (Model, ModelEntity, ModelEntity), g, ge1, ge2)

  return val != 0
end

"""
  Note that users must manually call `destroy` on the resulting model (ie.
  no finalizer is attached).
"""
function load(fname::String)

  g = ccall( (:gmi_load, gmi_lib), Ptr{Void}, (Cstring,), fname)
  return Model(g)
end

function destroy(g::Model)

  ccall( (:gmi_destroy, gmi_lib), Void, (Model,), g)
end

#------------------------------------------------------------------------------
# gmi_sim

function sim_start()
  ccall( (:gmi_sim_startJ, pumi_libname), Void, (), )
end

function sim_stop()
  ccall( (:gmi_sim_stopJ, pumi_libname), Void, (), )
end

function register_sim()
  ccall( (:gmi_register_simJ, pumi_libname), Void, (), )
end

# the functions gmi_sim_startJ/stopJ check if Simmetrix is supported, so call
# these unconditionally
sim_start()
#atexit(sim_stop)  # atexit hooks run before finalizers, so if a finalizer
#                  # causes the mesh to be destroyed, it won't be able to
#                  # destroy the geometric model (causes a segfault)

end # end module


text:
""""""
function heuristic(Fₕ, Fₗ=zeros(size(Fₕ)); transform=baselinetransform, filter=baselinefilter)
    transform = transform(Fₕ, Fₗ)
    Fₗ = Fₗ |> transform
    Fₕ = Fₕ |> transform
    filter = filter(Fₕ, Fₗ)
    Fₗ = Fₗ |> filter
    Fₕ = Fₕ |> filter
    scale = intervalscaled(Fₕ, Fₗ, filter=identity, transform=identity)
    Fₗ = Fₗ |> scale
    Fₕ = Fₕ |> scale
    function correction(F)
        Σₕ² = StatsBase.cov(Array(Fₕ), dims=2)
        𝐧 = sum(abs.(Array(Σₕ²)), dims=2)
        𝐧[𝐧 .== 0] .= Inf
        N⁻¹ = FeatureMatrix(inv(sqrt(Diagonal(𝐧[:]))), getnames(Fₕ))
        return FeatureMatrix(N⁻¹*F, getnames(F))
    end
    return F -> F |> transform |> filter |> scale |> correction
end
export heuristic


text:
###############################################################################
#
#   AutSymbol/ AutGroup / Automorphism
#
###############################################################################

struct RTransvect
    i::Int8
    j::Int8
end

struct LTransvect
    i::Int8
    j::Int8
end

struct FlipAut
    i::Int8
end

struct PermAut
    perm::Generic.perm{Int8}
end

struct Identity end

struct AutSymbol <: GSymbol
   id::Symbol
   pow::Int8
   fn::Union{LTransvect, RTransvect, PermAut, FlipAut, Identity}
end

mutable struct AutGroup{N} <: AbstractFPGroup
   objectGroup::FreeGroup
   gens::Vector{AutSymbol}
end

mutable struct Automorphism{N} <: GWord{AutSymbol}
    symbols::Vector{AutSymbol}
    modified::Bool
    savedhash::UInt
    parent::AutGroup{N}

    function Automorphism{N}(f::Vector{AutSymbol}) where {N}
        return new{N}(f, true, zero(UInt))
    end
end

export Automorphism, AutGroup, Aut, SAut

###############################################################################
#
#   Type and parent object methods
#
###############################################################################

elem_type(::AutGroup{N}) where N = Automorphism{N}

parent_type(::Automorphism{N}) where N = AutGroup{N}

###############################################################################
#
#   AutSymbol defining functions
#
###############################################################################

function (ϱ::RTransvect)(v, pow::Integer=1)
    @inbounds Groups.r_multiply!(v[ϱ.i], (v[ϱ.j]^pow).symbols, reduced=false)
    return v
end

function (λ::LTransvect)(v, pow::Integer=1)
    @inbounds Groups.l_multiply!(v[λ.i], (v[λ.j]^pow).symbols, reduced=false)
    return v
end

function (σ::PermAut)(v, pow::Integer=1)
   w = deepcopy(v)
   if pow == 1
       @inbounds for k in eachindex(v)
           v[k].symbols = w[σ.perm.d[k]].symbols
       end
   else
       s = (σ.perm^pow).d
       @inbounds for k in eachindex(v)
           v[k].symbols = w[s[k]].symbols
       end
   end
   return v
end

function (ɛ::FlipAut)(v, pow::Integer=1)
   @inbounds if isodd(pow)
       v[ɛ.i].symbols = inv(v[ɛ.i]).symbols
   end
   return v
end

(::Identity)(v, pow::Integer=1) = v

# taken from ValidatedNumerics, under under the MIT "Expat" License:
# https://github.com/JuliaIntervals/ValidatedNumerics.jl/blob/master/LICENSE.md
function subscriptify(n::Integer)
    subscript_0 = Int(0x2080) # Char(0x2080) -> subscript 0
    @assert 0 <= n <= 9
    return Char(subscript_0 + n)
    # return [Char(subscript_0 + i) for i in reverse(digits(n))])
end

function id_autsymbol()
   return AutSymbol(Symbol("(id)"), 0, Identity())
end

function rmul_autsymbol(i::Integer, j::Integer; pow::Integer=1)
    id = Symbol("ϱ", subscriptify(i), subscriptify(j))
    return AutSymbol(id, pow, RTransvect(i, j))
end

function lmul_autsymbol(i::Integer, j::Integer; pow::Integer=1)
    id = Symbol("λ", subscriptify(i), subscriptify(j))
    return AutSymbol(id, pow, LTransvect(i, j))
end

function flip_autsymbol(i::Integer; pow::Integer=1)
    if iseven(pow)
       return id_autsymbol()
    else
        id = Symbol("ɛ", subscriptify(i))
        return AutSymbol(id, 1, FlipAut(i))
    end
end

function perm_autsymbol(p::Generic.perm{I}; pow::Integer=one(I)) where I<:Integer
    if pow != 1
        p = p^pow
    end
    for i in eachindex(p.d)
        if p.d[i] != i
            id = Symbol("σ", [subscriptify(i) for i in p.d]...)
            return AutSymbol(id, 1, PermAut(p))
        end
    end
    return id_autsymbol()
end

function perm_autsymbol(a::Vector{T}) where T<:Integer
   return perm_autsymbol(perm(Vector{Int8}(a), false))
end

function domain(G::AutGroup{N}) where N
    F = G.objectGroup
    gg = gens(F)
    return ntuple(i->gg[i], Val(N))
end

###############################################################################
#
#   AutGroup / Automorphism constructors
#
###############################################################################

function AutGroup(G::FreeGroup; special=false)
   S = AutSymbol[]
   n = length(gens(G))
   n == 0 && return AutGroup{n}(G, S)

   indexing = [[i,j] for i in 1:n for j in 1:n if i≠j]

   rmuls = [rmul_autsymbol(i,j) for (i,j) in indexing]
   lmuls = [lmul_autsymbol(i,j) for (i,j) in indexing]

   append!(S, [rmuls; lmuls])

   if !special
      flips = [flip_autsymbol(i) for i in 1:n]
      syms = [perm_autsymbol(p) for p in PermutationGroup(n)][2:end]

      append!(S, [flips; syms])

   end
   return AutGroup{n}(G, S)
end

Aut(G::Group) = AutGroup(G)
SAut(G::Group) = AutGroup(G, special=true)

###############################################################################
#
#   Types call overloads
#
###############################################################################

Automorphism{N}(s::AutSymbol) where N = Automorphism{N}(AutSymbol[s])

function (G::AutGroup{N})() where N
   id = Automorphism{N}(id_autsymbol())
   id.parent = G
   return id
end

function (G::AutGroup{N})(f::AutSymbol) where N
   g = Automorphism{N}([f])
   g.parent = G
   return g
end

function (G::AutGroup{N})(g::Automorphism{N}) where N
   g.parent = G
   return g
end

###############################################################################
#
#   Functional call overloads for evaluation of AutSymbol and Automorphism
#
###############################################################################

function (f::AutSymbol)(v::NTuple{N, T}) where {N, T}
   if f.pow != 0
      v = f.fn(v, f.pow)::NTuple{N, T}
   end
   return v
end

function (F::Automorphism{N})(v::NTuple{N, T}) where {N, T}
    for f in reverse(F.symbols)
        v = f(v)::NTuple{N, T}
    end
    return v
end

###############################################################################
#
#   Comparison
#
###############################################################################

const HASHINGCONST = 0x7d28276b01874b19 # hash(Automorphism)

hash(s::AutSymbol, h::UInt) = hash(s.id, hash(s.pow, hash(:AutSymbol, h)))

function hash(g::Automorphism, h::UInt)
    if g.modified
        g_im = reduce!.(g(domain(parent(g))))
        g.savedhash = hash(g_im, hash(typeof(g), hash(parent(g), HASHINGCONST)))
        g.modified = false
        end
    return xor(g.savedhash, h)
end

function (==)(g::Automorphism{N}, h::Automorphism{N}) where N
    parent(g) == parent(h) || return false

    if !g.modified && !h.modified
        if g.savedhash != h.savedhash
            return false
        end
    end

    # expensive:
    g_im = reduce!.(g(domain(parent(g))))
    h_im = reduce!.(h(domain(parent(h))))
    # cheap:
    g.savedhash = hash(g_im, hash(typeof(g), hash(parent(g), HASHINGCONST)))
    g.modified = false
    h.savedhash = hash(h_im, hash(typeof(h), hash(parent(h), HASHINGCONST)))
    h.modified = false

    return g_im == h_im
end

###############################################################################
#
#   Basic manipulation
#
###############################################################################

function change_pow(s::AutSymbol, n::Integer)
    if n == zero(n)
        return id_autsymbol()
    end
    symbol = s.fn
    if symbol isa FlipAut
        return flip_autsymbol(symbol.i, pow=n)
    elseif symbol isa PermAut
        return perm_autsymbol(symbol.perm, pow=n)
    elseif symbol isa RTransvect
        return rmul_autsymbol(symbol.i, symbol.j, pow=n)
    elseif symbol isa LTransvect
        return lmul_autsymbol(symbol.i, symbol.j, pow=n)
    elseif symbol isa Identity
        return s
    else
        warn("Changing power of an unknown type of symbol! $s")
        return AutSymbol(s.id, n, s.fn)
    end
end

length(s::AutSymbol) = abs(s.pow)

###############################################################################
#
#   String I/O
#
###############################################################################

function show(io::IO, G::AutGroup)
   print(io, "Automorphism Group of $(G.objectGroup)\n")
   print(io, "Generated by $(join(G.gens, ","))")
end

###############################################################################
#
#   Binary operators
#
###############################################################################

###############################################################################
#
#   Inversion
#
###############################################################################

inv(f::AutSymbol) = change_pow(f, -f.pow)

###############################################################################
#
#   Misc
#
###############################################################################

function getperm(s::AutSymbol)
    if s.pow != 1
        @warn("Power for perm_symbol should be never 0!")
        return s.fn.perm^s.pow
    else
        return s.fn.perm
    end
end

function simplifyperms!(W::Automorphism{N}) where N
    reduced = true
    to_delete = Int[]
    for i in 1:length(W.symbols)-1
        if W.symbols[i].pow == 0
            continue
        elseif W.symbols[i].fn isa PermAut && W.symbols[i+1].fn isa PermAut
            reduced = false
            c = W.symbols[i]
            n = W.symbols[i+1]
            W.symbols[i+1] = perm_autsymbol(getperm(c)*getperm(n))
            push!(to_delete, i)
        end
    end
    deleteat!(W.symbols, to_delete)
    deleteids!(W)
    return reduced
end

function reduce!(W::Automorphism)
    if length(W) == 0
        return W
    elseif length(W.symbols) == 1
        deleteids!(W)
    else
        reduced = false
        while !reduced
            reduced = simplifyperms!(W) && freereduce!(W)
        end
    end

    W.modified = true

    return W
end

function linear_repr(A::Automorphism{N}, hom=matrix_repr) where N
    return reduce(*, linear_repr.(A.symbols, N, hom), init=hom(Identity(),N,1))
end

linear_repr(a::AutSymbol, n::Int, hom) = hom(a.fn, n, a.pow)

function matrix_repr(a::Union{RTransvect, LTransvect}, n::Int, pow)
    x = Matrix{Int}(I, n, n)
    x[a.i,a.j] = pow
    return x
end

function matrix_repr(a::FlipAut, n::Int, pow)
    x = Matrix{Int}(I, n, n)
    x[a.i,a.i] = -1^pow
    return x
end

matrix_repr(a::PermAut, n::Int, pow) = Matrix{Int}(I, n, n)[(a.perm^pow).d, :]

matrix_repr(a::Identity, n::Int, pow) = Matrix{Int}(I, n, n)


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

using Markdown
using InteractiveUtils

# ╔═╡ 5a0472c2-f40b-11ea-231b-a52e7382c1bc
md"""Recall the mean value theorem:

Let $f(x)$ be a function that is continuous on $[a,b]$ and differentiable on $(a,b)$.

Then there is a number $c$ in between $a,b$ such that:

$$f'(c) = \frac{f(b)-f(a)}{b-a}$$.

This means that there is a point with a derivative equal to the slope between the boundary points.

Let's take the following function:

$$f(x) = x^3+2x-x \text{ on } [-1,2]$$"""

# ╔═╡ bc85d0e0-f40d-11ea-2344-dd5ea4b66cde
md""" By the mean value theorem:

$$\frac{f(b)-f(a)}{b-a}= \frac{(10)-(-2)}{2-(-1)} = 4 = f'(c)$$

for some $c$ in $(-1,2)$.

With some algebra and rearranging, we find

$$c = \frac{-4+\sqrt{76}}{6}= 0.7863$$.

Graph:
"""

# ╔═╡ 40d54e00-f40b-11ea-3ddf-a39ce834ea41
using Gadfly

# ╔═╡ 54dca870-f411-11ea-1f67-ad216b5b2e49
set_default_plot_size(20cm, 15cm)

# ╔═╡ dd853160-f40c-11ea-318e-5751cdc6ffdf
begin
	p = plot((x)->(x^3+2x-x),-1,2)
	push!(p,Coord.cartesian(fixed=true))
	push!(p,Guide.title("f(x)"))
end

# ╔═╡ 5e773960-f40f-11ea-1894-0b4feb0e99e3
begin
	push!(p,layer(x=[0.7863],y=[1.2724], size=[0.25], color=[colorant"gold"]))
	push!(p,layer((x)->(4x-1.8728),0,2))
	push!(p,layer((x)->(4x+2),-1,2))
end

# ╔═╡ Cell order:
# ╟─40d54e00-f40b-11ea-3ddf-a39ce834ea41
# ╠═5a0472c2-f40b-11ea-231b-a52e7382c1bc
# ╟─54dca870-f411-11ea-1f67-ad216b5b2e49
# ╟─dd853160-f40c-11ea-318e-5751cdc6ffdf
# ╠═bc85d0e0-f40d-11ea-2344-dd5ea4b66cde
# ╠═5e773960-f40f-11ea-1894-0b4feb0e99e3


text:
#runs a regression on every firm
function regressbynonprofit!(df::DataFrame;
  YField::Symbol = :lreturn,
  bmfield::Symbol = Symbol(:pred_, YField),
  XFields::Vector{Symbol} = [:lsp500ret],
  groupbyfield::Symbol = :ein,
  suppressintercept::Bool = false,
  XNames::Vector{Symbol} = [:intercept; XFields],
  YName::Symbol = YField,
  βNames::Vector{Symbol} = (s::Symbol->Symbol(bmfield, "b", s)).(XNames),
  minpointspluskforβ::Int = MIN_POINTS_PLUS_K_FOR_BETA,
  predictwithresults::Bool = true,
  regstring::String = string(suppressintercept ? "0 + " : "", join((string).(XFields), " + ")),
  checkunique::Bool = false,
  winsorizeβ::Bool = false,
  winsorizeprop::Float64 = WINSORIZE_PROP
  )::Nothing

  #make sure we have the right number of names
  local N::Int = size(df,1)

  @assert length(XNames) == length(βNames)
  @assert length(XFields) == length(βNames) - 1 + suppressintercept

  local XSpec::FMExpr = Meta.parse(regstring)

  #preallcoate the arrays for the coefficients and prediction
  #println(skipmissing(df[1:1000, :pcontributions]))
  for f::Symbol ∈ βNames
    df[f] = Vector{MFloat64}(undef, N)
  end

  predictwithresults && (df[bmfield] = Vector{MFloat64}(undef, N))
  #split/apply/combine
  ctr=0

  #cache the fields we will use to test if there are enough data points
  checkfields::Vector{Symbol} = [YField; XFields;]

  #loop over all unique eins
  for subdf ∈ groupby(df, groupbyfield)

    local subN::Int
    if checkunique
      ssubdf = view(subdf, completecases(subdf[checkfields]), checkfields)
      subN = size(unique(ssubdf[XFields]),1)
    else
      subN = sum(completecases(subdf[checkfields]))
    end

    #make sure we have enough points for the regression
    if subN ≥ minpointspluskforβ + length(XNames)
      local lm::FMLM

      try
        lm = FMLM(subdf, XSpec, YField, XNames = XNames, YName=YName)
      catch err
        println(subdf[[:ein; :fisyr; checkfields]])
        error("ERROR: $err")
      end

      #predict the resutls and write down the coefficients
      for i ∈ 1:length(βNames)
        subdf[βNames[i]] .=  lm.β[i]
      end
    end
  end

  #note: for the perisstence results we do not want to winsorize
  if winsorizeβ
    for f ∈ βNames
      winsorize!(df[f], prop=winsorizeprop)
    end
  end

  if predictwithresults
    df[bmfield] .= suppressintercept ? 0.0 : df[βNames[1]]
    for i ∈ 1:length(XNames)
      if XNames[i] ≠ :intercept
        df[bmfield] .+= df[βNames[i]] .* df[XNames[i]]
      end
    end
  end

  return nothing
end

#helper function signature which takes in the data object
regressbynonprofit!(data::NCCSData;
  YField::Symbol = :lreturn,
  bmfield::Symbol = Symbol(:pred_, YField),
  XFields::Vector{Symbol} = [:lsp500ret],
  groupbyfield::Symbol = :ein,
  suppressintercept::Bool = false,
  XNames::Vector{Symbol} = [:intercept; XFields],
  YName::Symbol = YField,
  βNames::Vector{Symbol} = (s::Symbol->Symbol(bmfield, "b", s)).(XNames),
  minpointspluskforβ::Int = MIN_POINTS_PLUS_K_FOR_BETA,
  winsorizeβ::Bool = false) = (
    regressbynonprofit!(data.df,
      bmfield=bmfield,
      YField = YField,
      XFields = XFields,
      groupbyfield = groupbyfield,
      suppressintercept = suppressintercept,
      XNames = XNames,
      YName = YName,
      βNames = βNames,
      minpointspluskforβ = minpointspluskforβ,
      winsorizeβ=winsorizeβ))

function βdf(data::NCCSData; minpointspluskforβ::Int = MIN_POINTS_PLUS_K_FOR_BETA,
  path::String=WORKING_PATH, charitytypes::Symbol = :all)
  local βdf::DataFrame = DataFrame(ein = unique(data.df[:ein]))
  local Nβ::Int = size(βdf,1)
  local reg::FMLM

  local plotNames::Vector{String} = Vector{String}()
  local plots::Vector{PlotContainer} = Vector{PlotContainer}()

  #NOTE: DELETE THE BELOW

  (charitytypes ≠ :all) && (βdf = βdf[βdf[:charitytype] .== charitytypes, :])

  #initialize the rest of the dataframe
  βdf[:b0return] = Vector{MFloat64}(missing,Nβ)
  βdf[:b1return] = similar(βdf[:b0return])

  βdf[:b0contributions] = similar(βdf[:b0return])
  βdf[:b1contributions] = similar(βdf[:b0return])

  βdf[:b0contribonret] = similar(βdf[:b0return])
  βdf[:b1contribonret] = similar(βdf[:b0return])

  βdf[:b0expenses] = similar(βdf[:b0return])
  βdf[:b1expenses] = similar(βdf[:b0return])

  βdf[:avgwealth] = similar(βdf[:b0return])
  #βdf[:N] = Vector{MInt}(missing, Nβ)
  βdf[:charitytype] = Vector{MSymbol}(missing, Nβ)

  local YFields::Vector{Symbol} = Vector{Symbol}() # holds the LHS
  local XFields::Vector{Vector{Symbol}} = Vector{Vector{Symbol}}() #holds the RHS
  local βNames::Vector{Vector{Symbol}} = Vector{Vector{Symbol}}() # holds the names of the x variables

  #regress excess returns on excess s&p
  push!(XFields, [:sp500ret])
  push!(βNames, [:b0return, :b1return])
  #push!(XNames, [:intercept; :sp500ret])
  push!(YFields, :excess)

  #regress contributions on excess s&p
  push!(XFields, [:sp500ret])
  push!(βNames, [:b0contributions, :b1contributions])
  #push!(XNames, [:intercept; :lagsp500ret])
  push!(YFields, :pcontributions)

  #regress contributions on excess returns
  push!(XFields, [:excess])
  push!(βNames, [:b0contribonret, :b1contribonret])
  #push!(XNames, [:intercept; :excess])
  push!(YFields, :pcontributions)

  #regress contributions on excess returns
  push!(XFields, [:sp500ret])
  push!(βNames, [:b0expenses, :b1expenses])
  #push!(XNames, [:intercept; :sp500ret])
  push!(YFields, :pexpenses)

  println("beginning beta computations")
  for i ∈ 1:length(YFields)

    regressbynonprofit!(data.df, YField=YFields[i], XFields=XFields[i],
      βNames=βNames[i], minpointspluskforβ=minpointspluskforβ)
  end

  allβ::Vector{Symbol} = [βNames...;]

  #println(data.df[100:110,[:ein; allβ;]])
  #println("Complete cases of b0return: $(sum(!ismissing(data.df[:b0return])))")

  #βdf = join(βdf, data.df[[:ein; :charitytype; allβ;]], on = :ein, kind=:left)
  #index the dataframe
  βdfindex::Dict = Dict(βdf[i,:ein] => βdf[i,:] for i::Int  ∈ 1:Nβ)
  for subdf ∈ groupby(data.df, :ein) #compute the average wealth for each firm
    subN::Int = sum(completecases(subdf[[:excess, :sp500ret, :pcontributions, :pexpenses]]))
    βrow = βdfindex[subdf[1,:ein]]
    βrow[:charitytype] = subdf[1,:charitytype]

    if subN ≥ minpointspluskforβ
      for f ∈ allβ
        βrow[f] = subdf[1,f]
      end

      βrow[:avgwealth] =
        mean(skipmissing((ismissing).(subdf[:excess]) .* subdf[:adjnetassets]))
    end

  end

  #histogram for betas of returns
  βdf = βdf[completecases(βdf[[:b1return, :b0return]]),:]
  push!(plotNames, "return_beta_hist_$charitytypes")
  push!(plots, plot(βdf[(f->abs(f)<1.1).(βdf[:b1return]),:], x=:b1return,
    Guide.title("histogram of return betas"), Geom.histogram(bincount=25)))

  #histogram for betas of the contributions
  push!(plotNames, "return_beta_contrib_$charitytypes")
  βdf = βdf[completecases(βdf[[:b1contributions, :b0contributions]]),:]
  push!(plots, plot(βdf[(f->abs(f)<1.1).(βdf[:b1contributions]),:], x=:b1contributions,
    Guide.title("histogram of contribution betas"), Geom.histogram(bincount=25)))

  #histogram for betas of the contribonret
  βdf = βdf[completecases(βdf[[:b1contribonret, :b0contribonret]]),:]
  push!(plotNames, "return_beta_ret_on_contrib_$charitytypes")
  push!(plots,  plot(βdf[(f->abs(f)<1.1).(βdf[:b1contribonret]),:], x=:b1contribonret,
    Guide.title("histogram of return on contribution betas"), Geom.histogram(bincount=25)))

    #histogram for betas of expenses
  βdf = βdf[completecases(βdf[[:b1expenses, :b0contribonret]]),:]
  push!(plotNames, "return_beta_expense_$charitytypes")
  push!(plots, plot(βdf[(f->abs(f)<1.1).(βdf[:b1expenses]),:], x=:b1expenses,
      Guide.title("histogram of expense betas"), Geom.histogram(bincount=25)))

  println("return beta computations complete with median $(median(βdf[:b1return])),
    mean $(mean(βdf[:b1return])),
    weighted mean $(mean(βdf[:b1return], Weights(Vector{Float64}(βdf[:avgwealth])))),
    sd $(std(βdf[:b1return])), and se $(std(βdf[:b1return])/size(βdf,1)^0.5)")

  println("contribution on return beta computations complete with median $(median(βdf[:b1contributions])),
    mean $(mean(βdf[:b1contributions])),
    weighted mean $(mean(βdf[:b1contributions], Weights(Vector{Float64}(βdf[:avgwealth])))),
    sd $(std(βdf[:b1contributions])), and se $(std(βdf[:b1contributions])/size(βdf,1)^0.5)")

  println("contribution on return beta computations complete with median $(median(βdf[:b1contribonret])),
    mean $(mean(βdf[:b1contribonret])),
    weighted mean $(mean(βdf[:b1contribonret], Weights(Vector{Float64}(βdf[:avgwealth])))),
    sd $(std(βdf[:b1contribonret])), and se $(std(βdf[:b1contribonret])/size(βdf,1)^0.5)")

  println("expense beta computations complete with median $(median(βdf[:b1expenses])),
    mean $(mean(βdf[:b1expenses])),
    weighted mean $(mean(βdf[:b1expenses], Weights(Vector{Float64}(βdf[:avgwealth])))),
    sd $(std(βdf[:b1expenses])), and se $(std(βdf[:b1expenses])/size(βdf,1)^0.5)")

  bronbc = FMLM(βdf, Meta.parse("b1return"), :b1contributions,
    XNames = [:b1return], YName=:b1contributions)
  println("beta contributions on beta returns is $(bronbc.β)")
  #NOTE: THinking of a table which has something like the above stats with a pooled beta stat too

  for i ∈ 1:length(plots) #write the graphs
    try
      draw(SVG("$path\\$(plotNames[i]).svg", 9inch, 7inch), plots[i])
    catch err
      @error "ERROR writing $(plotNames[i]). Error message $err"
      println("Size: $(size(subdf))")
    end
    #println("Graph $(plotNames[i]) written.")
  end

  return βdf
end

function processβ(data::NCCSData)::NCCSData
  βdf(data)

  return data
end

max_stars_repo_name: clintonTE/CCA

max_issues_repo_name: clintonTE/CCA

max_forks_repo_name: clintonTE/CCA


text:
using BosonSampling
using Test
using JLD

### scattering ###

m = 5
n = 3

interf = RandHaar(m)
i = Input{RandomGramMatrix}(first_modes(n,m))
o = FockDetection(first_modes(n,m))
ev = Event(i,o,interf)

compute_probability!(ev)

ev.proba_params

### bunching ###

m = 5
n = 3

interf = RandHaar(m)
ib = Input{Bosonic}(first_modes(n,m))
ipd = Input{RandomGramMatrix}(first_modes(n,m))
subset_modes = first_modes(n,m)

typeof(subset_modes)

pb = full_bunching_probability(interf, ib, subset_modes)
ppd = full_bunching_probability(interf, ipd, subset_modes)

@test pb/ppd > 1. # this doesn't HAVE TO pass but will pass in nearly all
# cases

### Classical sampling ###

m = 16
n = 3

input = Input{Distinguishable}(first_modes(n,m))
interf = RandHaar(m)
out = classical_sampler(input=input, interf=interf)

## Cliffords sampler ###

m = 16
n = 3
input = Input{Bosonic}(first_modes(n,m))
interf = RandHaar(m)
out = cliffords_sampler(input=input, interf=interf)

## Noisy sampling ###

m = 16
n = 3

x = 0.8 # distinguishability
η = 0.8 # reflectivity

input = Input{OneParameterInterpolation}(first_modes(n,m), x)
interf = RandHaar(m)

out = noisy_sampler(input=input, reflectivity=η, interf=interf)

### MIS sampling ###

n = 8
m = n^2
starting_state = zeros(Int, m)
input_state = first_modes_array(n,m)

U = copy(rand_haar(m))

# generate a collisionless state as a starting point
starting_state = iterate_until_collisionless(() -> random_occupancy(n,m))

known_pdf(state) = process_probability_distinguishable(U, input_state, state)
target_pdf(state) = process_probability(U, input_state, state)
known_sampler = () -> iterate_until_collisionless(() -> classical_sampler(U, m, n)) # gives a classical sampler

# samples = metropolis_sampler(;target_pdf = target_pdf, known_pdf = known_pdf, known_sampler = known_sampler, starting_state = starting_state, n_iter = 100)

### Noisy distribution ###

n = 3
m = 6
x = 0.8
η = 0.8

input = Input{OneParameterInterpolation}(first_modes(n,m), x)
interf = RandHaar(m)

output_statistics = noisy_distribution(input=input, reflectivity=η, interf=interf)
p_exact = output_statistics[1]
p_approx = output_statistics[2]
p_sampled = output_statistics[3]

### Theoretical distribution ###

n = 3
m = 6

input = Input{Bosonic}(first_modes(n,m))
interf = RandHaar(m)
output_distribution = theoretical_distribution(input=input, interf=interf)

### Usage Interferometer ###

B = BeamSplitter(1/sqrt(2))
n = 2 # photon number
m = 2 # mode number
proba_bunching = Vector{Float64}(undef, 0)

for x = 0.0:0.01:1.0
    local input = Input{OneParameterInterpolation}(first_modes(n,m), x)
    p_theo = theoretical_distribution(input=input, interf=B)

    push!(proba_bunching, p_theo[2]) # store the probabilty to observe one photon in each mode
end

### subsets ###

s1 = Subset([1,1,0,0,0])
s2 = Subset([0,0,1,1,0])
s3 = Subset([1,0,1,0,0])

"subsets are not allowed to overlap"

# check_subset_overlap([s1,s2,s3]) will fail

### HOM tests: one mode ###

input_state = Input{Bosonic}(first_modes(n,m))

set1 = [1,0]
physical_interferometer = Fourier(m)
part = Partition([Subset(set1)])

(physical_indexes, pdf) = compute_probabilities_partition(physical_interferometer, part, input_state)

### HOM tests: mode1, mode2 ###

n = 2
m = 2

input_state = Input{Bosonic}(first_modes(n,m))

set1 = [1,0]
set2 = [0,1]
physical_interferometer = Fourier(m)
part = Partition([Subset(set1), Subset(set2)])

(physical_indexes, pdf) = compute_probabilities_partition(physical_interferometer, part, input_state)

print_pdfs(physical_indexes, pdf,n; partition_spans_all_modes = true, physical_events_only = true)

# for a single count

part_occ = PartitionOccupancy(ModeOccupation([1,1]),2,part)

compute_probability_partition_occupancy(physical_interferometer, part_occ, input_state)

# the same using an event

PartitionCount(part_occ)

o = PartitionCount(part_occ)
ev = Event(input_state, o, physical_interferometer)
############ need to change the constructor of Event

get_parametric_type(input_state)
get_parametric_type(o)

length(collect(typeof(o).parameters))
typeof(o)
collect(typeof(input_state).parameters)

### multiset for a random interferometer ###

m = 4
n = 3

inp = Input{Bosonic}(first_modes(n,m))
set1 = zeros(Int,m)
set2 = zeros(Int,m)
set1[1:2] .= 1
set2[3:4] .= 1

physical_interferometer = RandHaar(m)
part = Partition([Subset(set1), Subset(set2)])

(physical_indexes,  pdf) = compute_probabilities_partition(physical_interferometer, part, inp)
fourier_indexes = copy(physical_indexes)

print_pdfs(physical_indexes, pdf, n; physical_events_only = true, partition_spans_all_modes = true)
#print_pdfs(physical_indexes,  probas_fourier, n)

### partitions, subsets ###

n = 2
m = 5

s1 = Subset([1,1,0,0,0])
s2 = Subset([0,0,1,1,0])

part = Partition([s1,s2])
part_occ = PartitionOccupancy(ModeOccupation([2,0]),n,part)

i = Input{Bosonic}(first_modes(n,m))
o = PartitionCount(part_occ)
interf = RandHaar(m)
ev = Event(i,o,interf)

compute_probability!(ev)

### multiple counts probabilities ###

m = 10
n = 3
set1 = zeros(Int,m)
set2 = zeros(Int,m)
set1[1:2] .= 1
set2[3:4] .= 1

interf = RandHaar(m)
part = Partition([Subset(set1), Subset(set2)])

i = Input{Bosonic}(first_modes(n,m))
o = PartitionCountsAll(part)
ev = Event(i,o,interf)

compute_probability!(ev)


text:
# ---
# title: 386. Lexicographical Numbers
# id: problem386
# author: Indigo
# date: 2021-06-15
# difficulty: Medium
# categories:
# link: <https://leetcode.com/problems/lexicographical-numbers/description/>
# hidden: true
# ---
#
# Given an integer _n_ , return 1 - _n_ in lexicographical order.
#
# For example, given 13, return: [1,10,11,12,13,2,3,4,5,6,7,8,9].
#
# Please optimize your algorithm to use less time and space. The input size may
# be as large as 5,000,000.
#
#
## @lc code=start
using LeetCode

function lexical_order(n)
    res = Vector{Int}(undef, n)
    cur = 1
    for i in 1:n
        res[i] = cur
        if 10cur <= n
            cur *= 10
        else
            (cur >= n) && (cur ÷= 10)
            cur += 1
            while cur % 10 == 0
                cur ÷= 10
            end
        end
    end
    return res
end
## @lc code=end


text:
function part1(input)
    good_passwords = 0
    for line in readlines(input)
        n1, n2, c, password = match(r"(\d+)-(\d+) (\w): (\w+)", line).captures
        good_passwords += parse(Int, n1) <= count(==(c), split(password, "")) <= parse(Int, n2)
    end
    return good_passwords
end

function part2(input)
    good_passwords = 0
    for line in readlines(input)
        n1, n2, c, password = match(r"(\d+)-(\d+) (\w): (\w+)", line).captures
        m1, m2 = parse.(Int, (n1, n2))
        good_passwords += xor(password[m1:m1] == c, password[m2:m2] == c)
    end
    return good_passwords
end


text:
# By default, Julia/LLVM does not use fused multiply-add operations (FMAs).
# Since these FMAs can increase the performance of many numerical algorithms,
# we need to opt-in explicitly.
# See https://ranocha.de/blog/Optimizing_EC_Trixi for further details.
@muladd begin

function allocate_coefficients(mesh::AbstractMesh, equations, dg::DG, cache)
  # We must allocate a `Vector` in order to be able to `resize!` it (AMR).
  # cf. wrap_array
  zeros(eltype(cache.elements), nvariables(equations) * nnodes(dg)^ndims(mesh) * nelements(dg, cache))
end

@inline function wrap_array(u_ode::AbstractVector, mesh::AbstractMesh, equations, dg::DGSEM, cache)
  @boundscheck begin
    @assert length(u_ode) == nvariables(equations) * nnodes(dg)^ndims(mesh) * nelements(dg, cache)
  end
  # We would like to use
  #     reshape(u_ode, (nvariables(equations), ntuple(_ -> nnodes(dg), ndims(mesh))..., nelements(dg, cache)))
  # but that results in
  #     ERROR: LoadError: cannot resize array with shared data
  # when we resize! `u_ode` during AMR.
  #
  # !!! danger "Segfaults"
  #     Remember to `GC.@preserve` temporaries such as copies of `u_ode`
  #     and other stuff that is only used indirectly via `wrap_array` afterwards!

  # Currently, there are problems when AD is used with `PtrArray`s in broadcasts
  # since LoopVectorization does not support `ForwardDiff.Dual`s. Hence, we use
  # optimized `PtrArray`s whenever possible and fall back to plain `Array`s
  # otherwise.
  if LoopVectorization.check_args(u_ode)
    # This version using `PtrArray`s from StrideArrays.jl is very fast and
    # does not result in allocations.
    #
    # !!! danger "Heisenbug"
    #     Do not use this code when `@threaded` uses `Threads.@threads`. There is
    #     a very strange Heisenbug that makes some parts very slow *sometimes*.
    #     In fact, everything can be fast and fine for many cases but some parts
    #     of the RHS evaluation can take *exactly* (!) five seconds randomly...
    #     Hence, this version should only be used when `@threaded` is based on
    #     `@batch` from Polyester.jl or something similar. Using Polyester.jl
    #     is probably the best option since everything will be handed over to
    #     Chris Elrod, one of the best performance software engineers for Julia.
    PtrArray(pointer(u_ode),
             (StaticInt(nvariables(equations)), ntuple(_ -> StaticInt(nnodes(dg)), ndims(mesh))..., nelements(dg, cache)))
            #  (nvariables(equations), ntuple(_ -> nnodes(dg), ndims(mesh))..., nelements(dg, cache)))
  else
    # The following version is reasonably fast and allows us to `resize!(u_ode, ...)`.
    unsafe_wrap(Array{eltype(u_ode), ndims(mesh)+2}, pointer(u_ode),
                (nvariables(equations), ntuple(_ -> nnodes(dg), ndims(mesh))..., nelements(dg, cache)))
  end
end

# General fallback
@inline function wrap_array(u_ode::AbstractVector, mesh::AbstractMesh, equations, dg::DG, cache)
  wrap_array_native(u_ode, mesh, equations, dg, cache)
end

# Like `wrap_array`, but guarantees to return a plain `Array`, which can be better
# for interfacing with external C libraries (MPI, HDF5, visualization),
# writing solution files etc.
@inline function wrap_array_native(u_ode::AbstractVector, mesh::AbstractMesh, equations, dg::DG, cache)
  @boundscheck begin
    @assert length(u_ode) == nvariables(equations) * nnodes(dg)^ndims(mesh) * nelements(dg, cache)
  end
  unsafe_wrap(Array{eltype(u_ode), ndims(mesh)+2}, pointer(u_ode),
              (nvariables(equations), ntuple(_ -> nnodes(dg), ndims(mesh))..., nelements(dg, cache)))
end

function compute_coefficients!(u, func, t, mesh::AbstractMesh{1}, equations, dg::DG, cache)

  @threaded for element in eachelement(dg, cache)
    for i in eachnode(dg)
      x_node = get_node_coords(cache.elements.node_coordinates, equations, dg, i, element)
      u_node = func(x_node, t, equations)
      set_node_vars!(u, u_node, equations, dg, i, element)
    end
  end
end

function compute_coefficients!(u, func, t, mesh::AbstractMesh{2}, equations, dg::DG, cache)

  @threaded for element in eachelement(dg, cache)
    for j in eachnode(dg), i in eachnode(dg)
      x_node = get_node_coords(cache.elements.node_coordinates, equations, dg, i, j, element)
      u_node = func(x_node, t, equations)
      set_node_vars!(u, u_node, equations, dg, i, j, element)
    end
  end
end

function compute_coefficients!(u, func, t, mesh::AbstractMesh{3}, equations, dg::DG, cache)

  @threaded for element in eachelement(dg, cache)
    for k in eachnode(dg), j in eachnode(dg), i in eachnode(dg)
      x_node = get_node_coords(cache.elements.node_coordinates, equations, dg, i, j, k, element)
      u_node = func(x_node, t, equations)
      set_node_vars!(u, u_node, equations, dg, i, j, k, element)
    end
  end
end

end # @muladd


text:
# Searches for the (sub-)list containing a given needle within a provided list.
# If no such container can be found, an empty list is returned.
function container(haystack::Vector{Any}, needle::Any)
  # Search all sub-containers for the needle, recursively.
  cnt = {}
  for sub in haystack
    cnt = container(sub, needle)
    if !isempty(cnt)
      return cnt
    end
  end

  # Once we've finished searching the bottom of the haystack,
  # look at the container we're in.
  if in(needle, haystack)
    return haystack
  end

  # If we can't find the needle, then return the empty list.
  return {}
end

# If the haystack isn't a vector, then return an empty list.
container(haystack::Any, needle::Any) = {}

# Calculates the number of points in a given Push expression.
function num_points(expr::Vector{Any})
  points = 1
  for sub in expr
    points += num_points(sub)
  end
  return points
end
num_points(expr::Any) = 1

function insert_at_point!(v::Any, r::Any, pt::Integer)
  pt = abs(pt) % num_points(v)
  if pt == 0
    return r
  else
    insert_at_point!(v, r, pt, 0)
    return v
  end
end

function insert_at_point!(v::Vector{Any}, r::Any, pt::Integer, i::Integer)
  j = 0
  l = length(v)
  while i < pt && j < l
    j += 1
    i += 1
    if i == pt
      v[j] = r
    elseif isa(v[j], Vector)
      i = insert_at_point!(v[j], r, pt, i)
    end
  end
  return i
end

function position(haystack::Vector{Any}, needle::Any)
end


text:
abstract AbstractBot

# generic bot functions
is_alive(b::AbstractBot) = b.alive
team(b::AbstractBot) = b.team_id
botid(b::AbstractBot) = b.bot_id
name(b::AbstractBot) = b.name
# these might fail if they are not known
position(b::AbstractBot) = b.pos
hitpoints(b::AbstractBot) = b.HP

################################################
# bot types:
# own bot: we know everything about our own bots
type WebSockOwnBot <: AbstractBot
  bot_id::Int
  name::ASCIIString
  team_id::Int
  alive::Bool
  pos::Position
  HP::Int
end

# enemy bots: position and hp are unknown
type WebSockEnemyBot <: AbstractBot
  bot_id::Int
  name::ASCIIString
  team_id::Int
  alive::Bool
end

position(b::WebSockEnemyBot) = error("position of enemy bots are not known!")
hitpoints(b::WebSockEnemyBot) = error("hitpoints of enemy bots are not known!")

###################################################
# make a bot from the dict of the websocket message
function make_bot(d::Dict)
  # common data
  id = Int(d["botId"])
  name = d["name"]
  team = Int(d["teamId"])
  alive = Bool(d["alive"])

  if haskey(d, "pos")
    return WebSockOwnBot(id, name, team, alive, Position(d["pos"]), d["hp"])
  else
    return WebSockEnemyBot(id, name, team, alive)
  end
end

max_stars_repo_name: Condimentum/spacetyckiting

max_issues_repo_name: Condimentum/spacetyckiting

max_forks_repo_name: Condimentum/spacetyckiting


text:
using Documenter, ModiaPlot

makedocs(
  #modules  = [ModiaPlot],
  sitename = "ModiaPlot",
  authors  = "Martin Otter (DLR-SR)",
  format = Documenter.HTML(prettyurls = false),
  pages    = [
     "Home"      => "index.md",
	 "Functions" => "Functions.md",
  	 "Internal"  => "Internal.md",
  ]
)