text:
##############################################################################
##
## sumby...
##
##############################################################################

"""
Perform sum by group
```julia
sumby(df::Union{AbstractDataFrame,NDSparse}, by::Symbol, val::Symbol)
sumby(by::AbstractVector  val::AbstractVector)
```
### Arguments
* `df` : an AbstractDataFrame/NDSparse from which to extract the by and val columns
* `by` : data table column to group by
* `val`: data table column to sum
### Returns
* `::Dict` : A Dict that maps unqiues values of by to sum of val

### Examples
```julia
using FastGroupBy
using DataFrames, IndexedTables, Compat, BenchmarkTools
import DataFrames.DataFrame

const N = 10_000_000; const K = 100

# sumby is faster for DataFrame without missings
srand(1);
@time df = DataFrame(id = rand(1:Int(round(N/K)), N), val = rand(round.(rand(K)*100,4), N));
@belapsed DataFrames.aggregate(df, :id, sum)
@belapsed sumby(df, :id, :val)

```
"""
sumby(by, val, alg = :auto) = sumby!(copy(by), copy(val), alg)

# sumby(by, val) = sumby!(copy(by), copy(val))

function sumby!(by::AbstractVector{T},  val::AbstractVector{S}; alg = :auto)::Dict{T,S} where {T, S<:Number}
    l = length(by)

    l == length(val) || throw(ErrorException("length of by and val must be the same"))

    if l == 0
        return Dict{T,S}()
    elseif l == 1
        return Dict{T,S}(by[1], val[1])
    elseif issorted(by)
        return sumby_contiguous(by, val)
    elseif l <= 2^16
        return sumby_sortperm(by, val)
    elseif !isbits(T) || alg == :dict
        return sumby_dict(by, val)
    elseif nthreads() > 1
        return sumby_multi_rs(by, val)
    elseif l <= 50_000_000
        return sumby_radixsort!(by, val)
    else
        # return sumby_radixgroup(by, val)
        return sumby_radixsort!(by, val)
    end
end

"sumby by using radix and counting sort to group by; it's only a partial sort. It's faster for large by"
function sumby_radixgroup!(by::AbstractVector{T},  val::AbstractVector{S}; cutsize = 2048)::Dict{T,S} where {T, S<:Number}
    by_sim = similar(by)
    val1=similar(val)
    o = Forward

    l = length(by)

    if l == 0
        return Dict{T,S}()
    elseif l == 1
        return Dict{T,S}(by[1] => val[1])
    end

    # Make sure we're sorting a bits type
    if !isbits(T)
      error("sumby_nosort on works on bits types (got $T)")
    end

    # Init
    iters = ceil(Integer, sizeof(T)*8/RADIX_SIZE)

    # distinct range to groupby over; initially the whole block is one range
    nextranges = Tuple{Int, Int}[]

    new_small_grps = similar(nextranges)

    push!(nextranges, (1,l))

    for j = 1:iters
        #println(string(j,":",length(nextranges)))
        ranges = nextranges::Array{Tuple{Int, Int}}
        nextranges = Tuple{Int,Int}[]
        if length(ranges) > 0
            for (lo, hi) in ranges
                bin = zeros(UInt32, 2^RADIX_SIZE)

                # Histogram for each element, radix
                for i = lo:hi
                    #v = uint_mapping(o, by[i])
                    idx = Int((by[i] >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1
                    @inbounds bin[idx] += 1
                end

                # Sort!
                len = hi - lo + 1

                # Unroll first data iteration, check for degenerate case
                #v = uint_mapping(o, by[hi])
                idx = Int((by[hi] >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1

                # are all values the same at this radix?
                if bin[idx] == len
                    nextranges = vcat(nextranges, [(lo, hi)])::Array{Tuple{Int, Int}}
                    @inbounds by_sim[lo:hi] = by[lo:hi]
                    @inbounds val1[lo:hi] = val[lo:hi]
                else
                    if lo > 1;  bin[1] += UInt32(lo-1);  end
                    cbin = cumsum(bin)

                    # compute the new ranges for next round
                    tmpr = sort(setdiff(unique(cbin), lo-1))::Vector{Int}
                    tmpnextranges = vcat(nextranges, [(lo,Int(tmpr[1]))], [(Int(a)+1,Int(b)) for (a,b) in zip(tmpr[1:end-1],tmpr[2:end])])::Array{Tuple{Int, Int}}
                    nextranges = filter(x -> x[2] - x[1] + 1 > cutsize, tmpnextranges)::Array{Tuple{Int, Int}}

                    if length(tmpnextranges) != length(nextranges)
                        # this is a small grp just sort it
                        new_small_grps = filter(x -> x[2] - x[1] + 1 <= cutsize, tmpnextranges)
                    end

                    ci = cbin[idx]
                    by_sim[ci] = by[hi]
                    val1[ci] = val[hi]
                    cbin[idx] -= 1

                    # Finish the loop...
                    @inbounds for i in hi-1:-1:lo
                      #v = uint_mapping(o, by[i])
                      idx =Int((by[i] >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1
                      ci = cbin[idx]
                      by_sim[ci] = by[i]
                      val1[ci] = val[i]
                      cbin[idx] -= 1
                    end
                end
            end
            if length(new_small_grps) > 0
                for (lo_nsg, hi_nsg) in new_small_grps
                    @inbounds bytmp = by_sim[lo_nsg:hi_nsg]
                    sp = sortperm(bytmp)
                    bytmp_sorted = bytmp[sp]
                    val1_sorted = val1[lo_nsg:hi_nsg][sp]
                    @inbounds by[lo_nsg:hi_nsg] = bytmp_sorted
                    @inbounds val[lo_nsg:hi_nsg] = val1_sorted
                    @inbounds by_sim[lo_nsg:hi_nsg] = bytmp_sorted
                    @inbounds val1[lo_nsg:hi_nsg] = val1_sorted
                end
                new_small_grps = Int[]
            end
            by, by_sim = by_sim, by
            val, val1 = val1, val
        end
    end

    return sumby_contiguous(by, val)
end

"sumby by sorting the by column using radixsort"
function sumby_radixsort!(by::AbstractVector{T},  val::AbstractVector{S})::Dict{T,S} where {T, S<:Number}

  hi = length(by)

  if hi == 0;  return Dict{T,S}();
  elseif hi == 1;  return Dict{T,S}(by[1] => val[1]);  end

  by_sim = similar(by)
  val1=similar(val)
  lo = 1

  # Make sure we're sorting a bits type
  #TT = Base.Order.ordtype(o, by)
  if !isbits(T)
      throw(error("Radix sort only sorts bits types (got $T)"))
  end

  # Init
  iters = ceil(Integer, sizeof(T)*8/RADIX_SIZE)
  bin = zeros(UInt32, 2^RADIX_SIZE, iters)
  if lo > 1;  bin[1,:] = lo-1;  end

  # Histogram for each element, radix
  for i = lo:hi
    v = by[i]
    for j = 1:iters
        idx = Int((v >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1
        @inbounds bin[idx,j] += 1
    end
  end

  # Sort!
  # swaps = 0
  len = hi-lo+1
  for j = 1:iters
      # Unroll first data iteration, check for degenerate case
      v = by[hi]
      idx = Int((v >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1

      # are all values the same at this radix?
      if bin[idx,j] == len;  continue;  end

      cbin = cumsum(bin[:,j])
      ci = cbin[idx]
      by_sim[ci] = by[hi]
      val1[ci] = val[hi]

      cbin[idx] -= 1

      # Finish the loop...
      @inbounds for i in hi-1:-1:lo
          v = by[i]
          idx = Int((v >> (j-1)*RADIX_SIZE) & RADIX_MASK) + 1
          ci = cbin[idx]
          by_sim[ci] = by[i]
          val1[ci] = val[i]
          cbin[idx] -= 1
      end
      by,by_sim = by_sim,by
      val,val1 = val1,val
    #   swaps += 1
  end

  # @inbounds if isodd(swaps)
  #   by,by_sim = by_sim,by
  #   val,val1 = val1,val
  #   for i = lo:hi
  #     by[i] = by_sim[i]
  #     val[i] = val1[i]
  #   end
  # end

  sumby_contiguous(by, val)
end

"sumby assuming that the elements are organised contiguously; it does not perform a check"
function sumby_contiguous(by_sorted::AbstractVector{T},  val::AbstractVector{S})::Dict{T,S} where {T, S<:Number}
  res = Dict{T,S}()
  @inbounds tmp_val = val[1]
  @inbounds last_byi = by_sorted[1]
  @inbounds for i in 2:length(by_sorted)
    if by_sorted[i] == last_byi
      tmp_val += val[i]
    else
      res[last_byi] = tmp_val
      tmp_val = val[i]
      last_byi = by_sorted[i]
    end
  end

  @inbounds res[last_byi] = tmp_val

  res
end

#sumby_contiguous2 is too slow
# function sumby_contiguous2{T,S}(by_sorted::AbstractVector{T},  val::AbstractVector{S})
#   res = Dict{T,S}()
#   @inbounds last_byi = by_sorted[1]
#   lo = 1
#   hi = 1
#   @inbounds for i in 2:length(by_sorted)
#     if val[i] == last_byi
#       hi += 1
#     else
#       res[last_byi] = sum(val[lo:hi])
#       last_byi = by_sorted[i]
#       lo, hi = hi + 1, hi + 1
#     end
#   end
#
#   @inbounds res[last_byi] = sum(val[lo:hi])
#
#   res
# end

"This is faster for smaller by and also doesn't change the input"
function sumby_sortperm(by::AbstractVector{T},  val::AbstractVector{S})::Dict{T,S} where {T, S<:Number}
    sp = sortperm(by, alg = RadixSort)
    sumby_contiguous(view(by, sp), view(val,sp))
end
# function sumby_sortperm2{T, S<:Number}(by::AbstractVector{T},  val::AbstractVector{S})::Dict{T,S}
#     sp = sortperm(by)
#     sumby_contiguous(by[sp], val[sp])
# end

"sumby using Dict - can be quite slow due to slow hash table operations"
function sumby_dict(by::AbstractVector{T}, val::AbstractVector{S})::Dict{T,S} where {T,S<:Number}
  res = Dict{T, S}()
  # resize the Dict to a larger size
  for (byi, vali) in zip(by, val)
    index = ht_keyindex2!(res, byi)
    if index > 0
      @inbounds res.vals[index] += vali
    else
      @inbounds _setindex!(res, vali, byi, -index)
    end
  end
  return res
end

#Optimized sumby for PooledArrays
function sumby!(by::Union{PooledArray, CategoricalArray}, val::AbstractVector{S}) where {S<:Number}
  l = length(by.pool)
  res = zeros(S, l)

  for (i, v) in zip(by.refs, val)
    @inbounds res[i] += v
  end
  return Dict(by.pool[i] => res[i] for i in S(1):S(l))
end

sumby(by::Union{PooledArray, CategoricalArray}, val::AbstractVector{S}) where {S<:Number} = sumby!(by, val)

sumby!(dt::Union{AbstractDataFrame, NDSparse}, by::Symbol, val::Symbol) = sumby!(column(dt,by), column(dt,val))


text:
using Distributions

export EditAllocationMove

struct EditAllocationMove <: InvMcmcMove
    ν::Int
    counts::Vector{Int}  # Track acceptance
    function EditAllocationMove(;ν::Int=3)
        new(
            ν,
            # zeros(Int,ν),
            # [0],
            [0,0]
        )
    end
end

Base.show(io::IO, x::EditAllocationMove) = print(io, "EditAllocationMove(ν=$(x.ν))")

function prop_sample!(
    S_curr::InteractionSequence{Int},
    S_prop::InteractionSequence{Int},
    move::EditAllocationMove,
    pointers::InteractionSequence{Int},
    V::UnitRange
    )

    N = length(S_curr)
    δ = rand(1:move.ν)
    rem_edits = δ # Remaining edits to allocate
    len_diffs = 0
    # j = 0 # Keeps track how many interaction have been edited
    log_prod_term = 0.0

    for i in eachindex(S_curr)

        # If at end we just assign all remaining edits to final interaction
        if i == N
            δ_tmp = rem_edits
        # Otherwise we sample the number of edits via rescaled Binomial
        else
            p = 1/(N-i+1)
            δ_tmp = rand(Binomial(rem_edits, p)) # Number of edits to perform on ith interaction
        end

        if δ_tmp == 0
            continue
        else
            # j += 1 # Increment j
            # Make edits ....
            @inbounds n = length(S_curr[i])
            d = rand(0:min(n-1,δ_tmp))
            m = n + δ_tmp - 2*d

            I_tmp = S_prop[i]
            rand_delete!(I_tmp, d)
            rand_insert!(I_tmp, δ_tmp-d, V)

            # push!(mode.ind_upd, i)
            # @inbounds move.ind_upd[j] = i # Store which interaction was updated

            # Add to log_ratio
            # log_prod_term += log(b - a + 1) - log(ub(m, δ_tmp) - lb(m, δ_tmp, model) +1)
            log_prod_term += log(min(n-1, δ_tmp)+1) - log(min(m-1, δ_tmp)+1)
            len_diffs += m-n  # How much bigger the new interaction is
        end
        # move.num_upd[1] = j # Store number of paths updated

        # Update rem_edits
        rem_edits -= δ_tmp

        # If no more left terminate
        if rem_edits == 0
            break
        end

    end

    # # Add final part of log_ratio term
    log_ratio = log(length(V)) * len_diffs + log_prod_term

    return log_ratio

end

# function enact_accept!(
#     S_curr::InteractionSequence{Int},
#     S_prop::InteractionSequence{Int},
#     pointers::InteractionSequence{Int},
#     move::EditAllocationMove
#     )
#     ind = view(move.ind_upd, 1:move.num_upd[1])
#     for i in ind
#         @inbounds copy!(S_curr[i], S_prop[i])
#     end
# end

# function enact_reject!(
#     S_curr::InteractionSequence{Int},
#     S_prop::InteractionSequence{Int},
#     pointers::InteractionSequence{Int},
#     move::EditAllocationMove
#     )
#     ind = view(move.ind_upd, 1:move.num_upd[1])
#     for i in ind
#         @inbounds copy!(S_prop[i], S_curr[i])
#     end
# end

ext: jl


text:
#############################################################
#
# Favorites Functions
#
#############################################################

function get_favorites_list(; options=Dict{AbstractString, AbstractString}())

    r = get_oauth("https://api.twitter.com/1.1/favorites/list.json", options)

    #If 200 response, create Array{TWEETS,1}, otherwise return raw response
    return r.status == 200 ? to_TWEETS(Requests.json(r)) : r

end

function post_favorites_destroy(id::AbstractString; options=Dict{AbstractString, AbstractString}())

    options["id"] = id

    r = post_oauth("https://api.twitter.com/1.1/favorites/destroy.json", options)

    #If 200 response, create Array{TWEETS,1}, otherwise return raw response
    return r.status == 200 ? to_TWEETS(Requests.json(r)) : r

end

function post_favorites_create(id::AbstractString; options=Dict{AbstractString, AbstractString}())

    options["id"] = id

    r = post_oauth("https://api.twitter.com/1.1/favorites/create.json", options)

    #If 200 response, create Array{TWEETS,1}, otherwise return raw response
    return r.status == 200 ? to_TWEETS(Requests.json(r)) : r

end


text:
const LQOI = LinQuadOptInterface
const MOI  = LQOI.MOI

mutable struct LinQuadSOS
    kind::Symbol
    weights::Vector{Float64}
    indices::Vector{Int}
end

mutable struct MockLinQuadModel # <: LinQuadOptInterface.LinQuadOptimizer

    sense::Symbol

    Qcon::Vector{Matrix{Float64}}
    A::Matrix{Float64}
    b::Vector{Float64}
    c::Vector{Float64}
    Qobj::Matrix{Float64}

    range::Vector{Float64} # constraint UB for range

    lb::Vector{Float64}
    ub::Vector{Float64}

    vartype::Vector{Cchar} # :Bin, :Int, :Con
    contype::Vector{Cchar} # :LEQ :GEQ, :EQ, :RANGE

    sos::Vector{LinQuadSOS}

    termination_status::MOI.TerminationStatusCode
    termination_status_stored::Vector{MOI.TerminationStatusCode}

    primal_status::MOI.ResultStatusCode
    primal_status_stored::Vector{MOI.ResultStatusCode}

    dual_status::MOI.ResultStatusCode
    dual_status_stored::Vector{MOI.ResultStatusCode}

    variable_primal_solution::Vector{Float64}
    variable_primal_solution_stored::Vector{Vector{Float64}}
    variable_dual_solution::Vector{Float64}
    variable_dual_solution_stored::Vector{Vector{Float64}}

    constraint_primal_solution::Vector{Float64}
    constraint_primal_solution_stored::Vector{Vector{Float64}}
    constraint_dual_solution::Vector{Float64}
    constraint_dual_solution_stored::Vector{Vector{Float64}}

    quadratic_primal_solution::Vector{Float64}
    quadratic_primal_solution_stored::Vector{Vector{Float64}}

    quadratic_dual_solution::Vector{Float64}
    quadratic_dual_solution_stored::Vector{Vector{Float64}}

    ray_dual_solution::Vector{Float64}
    ray_dual_solution_stored::Vector{Vector{Float64}}

    ray_primal_solution::Vector{Float64}
    ray_primal_solution_stored::Vector{Vector{Float64}}

    function MockLinQuadModel(env::Void,str::String="defaultname")
        m = new()

        m.sense = :minimize

        m.Qcon = Matrix{Float64}[]
        m.A = zeros(0,0)
        m.b = zeros(0)
        m.c = zeros(0)
        m.Qobj = zeros(0,0)

        m.range = zeros(0)

        m.lb = zeros(0)
        m.ub = zeros(0)

        m.vartype = Cchar[]
        m.contype = Cchar[]

        m.sos = LinQuadSOS[]

        m.termination_status = MOI.Success
        m.termination_status_stored = MOI.TerminationStatusCode[]
        m.primal_status = MOI.FeasiblePoint
        m.primal_status_stored = MOI.ResultStatusCode[]
        m.dual_status = MOI.FeasiblePoint
        m.dual_status_stored = MOI.ResultStatusCode[]

        m.variable_primal_solution = zeros(0)
        m.variable_primal_solution_stored = Vector{Float64}[]
        m.variable_dual_solution = zeros(0)
        m.variable_dual_solution_stored = Vector{Float64}[]

        m.constraint_primal_solution = zeros(0)
        m.constraint_primal_solution_stored = Vector{Float64}[]
        m.constraint_dual_solution = zeros(0)
        m.constraint_dual_solution_stored = Vector{Float64}[]

        m.quadratic_primal_solution = zeros(0)
        m.quadratic_primal_solution_stored = Vector{Float64}[]
        m.quadratic_dual_solution = zeros(0)
        m.quadratic_dual_solution_stored = Vector{Float64}[]

        m.ray_primal_solution = zeros(0)
        m.ray_primal_solution_stored = Vector{Float64}[]
        m.ray_dual_solution = zeros(0)
        m.ray_dual_solution_stored = Vector{Float64}[]

        return m
    end
end

function unload(from,to,warn = true)
    if !isempty(from)
        out = from[1]
        shift!(from)
        return out
    else
        if warn
            warn("cant solve this model no extra solution")
        end
        return to
    end
end

function fakesolve(m::MockLinQuadModel)

    m.termination_status = unload(m.termination_status_stored,  m.termination_status)
    m.primal_status = unload(m.primal_status_stored,  m.primal_status)
    m.dual_status = unload(m.dual_status_stored,  m.dual_status)

    m.variable_primal_solution = unload(m.variable_primal_solution_stored,  m.variable_primal_solution)
    m.variable_dual_solution = unload(m.variable_dual_solution_stored,  m.variable_dual_solution)

    m.constraint_primal_solution = unload(m.constraint_primal_solution_stored,  m.constraint_primal_solution)
    m.constraint_dual_solution = unload(m.constraint_dual_solution_stored,  m.constraint_dual_solution)

    m.quadratic_primal_solution = unload(m.quadratic_primal_solution_stored,  m.quadratic_primal_solution)
    m.quadratic_dual_solution = unload(m.quadratic_dual_solution_stored,  m.quadratic_dual_solution)

    m.ray_primal_solution = unload(m.ray_primal_solution_stored,  m.ray_primal_solution, false)
    m.ray_dual_solution = unload(m.ray_dual_solution_stored,  m.ray_dual_solution, false)

    nothing
end

num_vars(inner::MockLinQuadModel) = length(inner.c)
num_cons(inner::MockLinQuadModel) = length(inner.b)

function set_variable_primal_solution!(inner::MockLinQuadModel,input)
    push!(m.variable_primal_solution_stored, input)
end
function set_variable_dual_solution!(inner::MockLinQuadModel,input)
    push!(m.variable_dual_solution_stored, input)
end
function set_constraint_primal_solution!(inner::MockLinQuadModel,input)
    push!(m.constraint_primal_solution_stored, input)
end
function set_constraint_dual_solution!(inner::MockLinQuadModel,input)
    push!(m.constraint_dual_solution_stored, input)
end
function set_quadratic_dual_solution!(inner::MockLinQuadModel,input)
    push!(m.quadratic_dual_solution_stored, input)
end
function set_quadratic_primal_solution!(inner::MockLinQuadModel,input)
    push!(m.quadratic_primal_solution_stored, input)
end
function set_ray_primal_solution!(inner::MockLinQuadModel,input)
    push!(m.ray_primal_solution_stored, input)
end
function set_ray_dual_solution!(inner::MockLinQuadModel,input)
    push!(m.ray_dual_solution_stored, input)
end
function set_termination_status!(inner::MockLinQuadModel,input)
    push!(m.termination_status_stored, input)
end
function set_primal_status!(inner::MockLinQuadModel,input)
    push!(m.primal_status_stored, input)
end
function set_dual_status!(inner::MockLinQuadModel,input)
    push!(m.dual_status_stored, input)
end

const SUPPORTED_OBJECTIVES = [
    LQOI.Linear,
    LQOI.Quad
]

const SUPPORTED_CONSTRAINTS = [
    (LQOI.Linear, LQOI.EQ),
    (LQOI.Linear, LQOI.LE),
    (LQOI.Linear, LQOI.GE),
    (Linear, IV),
    (LQOI.Quad, LQOI.EQ),
    (LQOI.Quad, LQOI.LE),
    (LQOI.Quad, LQOI.GE),
    (LQOI.SinVar, LQOI.EQ),
    (LQOI.SinVar, LQOI.LE),
    (LQOI.SinVar, LQOI.GE),
    (LQOI.SinVar, LQOI.IV),
    (LQOI.SinVar, MOI.ZeroOne),
    (LQOI.SinVar, MOI.Integer),
    (LQOI.VecVar, LQOI.SOS1),
    (LQOI.VecVar, LQOI.SOS2),
    (LQOI.VecVar, MOI.Nonnegatives),
    (LQOI.VecVar, MOI.Nonpositives),
    (LQOI.VecVar, MOI.Zeros),
    (LQOI.VecLin, MOI.Nonnegatives),
    (LQOI.VecLin, MOI.Nonpositives),
    (LQOI.VecLin, MOI.Zeros)
]

mutable struct MockLinQuadOptimizer <: LQOI.LinQuadOptimizer
    LQOI.@LinQuadOptimizerBase
    params::Dict{String,Any}
    l_rows::Vector{Int}
    q_rows::Vector{Int}
    MockLinQuadOptimizer(::Void) = new()
end

LQOI.LinearQuadraticModel(::Type{MockLinQuadOptimizer},env::Void) = MockLinQuadModel(env,"defaultname")

function MockLinQuadOptimizer(;kwargs...)

    m = MockLinQuadOptimizer(nothing)
    m.params = Dict{String,Any}()
    MOI.empty!(m)
    for (name,value) in kwargs
        m.params[string(name)] = value
        setparam!(m.inner, string(name), value)
    end
    return m
end

LQOI.get_number_linear_constraints(instance::MockLinQuadOptimizer) = num_cons(instance.inner) - LQOI.get_number_quadratic_constraints(instance)
function LQOI.get_number_quadratic_constraints(instance::MockLinQuadOptimizer)
    c = 0
    for i in 1:num_cons(instance.inner)
        if !isempty(instance.inner.Qcon[i])
            c += 1
        end
    end
    return c
end

function MOI.empty!(m::MockLinQuadOptimizer)
    MOI.empty!(m,nothing)
    m.l_rows = Int[]
    m.q_rows = Int[]
    for (name,value) in m.params
        setparam!(m.inner, name, value)
    end
end

LQOI.supported_constraints(s::MockLinQuadOptimizer) = SUPPORTED_CONSTRAINTS
LQOI.supported_objectives(s::MockLinQuadOptimizer) = SUPPORTED_OBJECTIVES

cintvec(v::Vector) = convert(Vector{Int32}, v)

backend_type(m::MockLinQuadOptimizer, ::MOI.GreaterThan{T}) where T = Cchar('G')
backend_type(m::MockLinQuadOptimizer, ::MOI.LessThan{T}) where T    = Cchar('L')
backend_type(m::MockLinQuadOptimizer, ::MOI.EqualTo{T}) where T     = Cchar('E')
# Implemented separately
# backend_type(m::MockLinQuadOptimizer, ::MOI.Interval{T}) where T    = Cchar('R')

backend_type(m::MockLinQuadOptimizer, ::MOI.Zeros)        = Cchar('E')
backend_type(m::MockLinQuadOptimizer, ::MOI.Nonpositives) = Cchar('L')
backend_type(m::MockLinQuadOptimizer, ::MOI.Nonnegatives) = Cchar('G')

backend_type(m::MockLinQuadOptimizer, ::Val{:Continuous}) = Cchar('C')
backend_type(m::MockLinQuadOptimizer, ::Val{:Upperbound}) = Cchar('U')
backend_type(m::MockLinQuadOptimizer, ::Val{:Lowerbound}) = Cchar('L')

# TODO - improve single type
function LQOI.change_variable_bounds!(instance::MockLinQuadOptimizer, colvec, valvec, sensevec)

    lb_len = count(x->x==Cchar('L'), sensevec)
    LB_val = Array{Float64}(0)
    sizehint!(LB_val, lb_len)
    LB_col = Array{Cint}(0)
    sizehint!(LB_col, lb_len)

    ub_len = count(x->x==Cchar('U'), sensevec)
    UB_val = Array{Float64}(0)
    sizehint!(UB_val, ub_len)
    UB_col = Array{Cint}(0)
    sizehint!(UB_col, ub_len)

    if lb_len + ub_len != length(sensevec)
        error("Invalid values for sensevec")
    end

    for i in eachindex(valvec)
        if sensevec[i] == Cchar('L')
            push!(LB_col, colvec[i])
            push!(LB_val, valvec[i])
        elseif sensevec[i] == Cchar('U')
            push!(UB_col, colvec[i])
            push!(UB_val, valvec[i])
        end
    end

    if lb_len > 0
        for i in eachindex(LB_col)
            instance.inner.lb[LB_col[i]] = LB_val[i]
        end
    end

    if ub_len > 0
        for i in eachindex(UB_col)
            instance.inner.ub[UB_col[i]] = UB_val[i]
        end
    end

    nothing
end

function LQOI.get_variable_lowerbound(instance::MockLinQuadOptimizer, col)
    instance.inner.lb[col]
end

function LQOI.get_variable_upperbound(instance::MockLinQuadOptimizer, col)
    instance.inner.ub[col]
end

function LQOI.add_linear_constraints!(instance::MockLinQuadOptimizer, A::CSRMatrix{Float64}, sensevec, rhsvec)

    newrows = length(rhsvec)
    rows = num_cons(instance.inner)
    addedrows = collect((rows+1):(rows+newrows))
    # quadind
    append!(instance.l_rows,addedrows)

    rowvec, colvec, coefvec = A.row_pointers, A.columns, A.coefficients

    rows = length(rhsvec)
    cols = size(instance.inner.A)[2]
    push!(rowvec,length(colvec)+1)

    An = full(SparseMatrixCSC(cols,rows,rowvec,colvec,coefvec)')

    instance.inner.A = vcat(instance.inner.A,An)

    append!(instance.inner.b,rhsvec)
    append!(instance.inner.contype,sensevec)

    append!(instance.inner.range,[Inf for i in 1:rows])
    append!(instance.inner.Qcon,[zeros(0,0) for i in 1:rows])

    nothing
end

function LQOI.add_ranged_constraints!(instance::MockLinQuadOptimizer, A::CSRMatrix{Float64}, lowerbound, upperbound)

    newrows = length(lowerbound)
    rows = num_cons(instance.inner)
    addedrows = collect((rows+1):(rows+newrows))
    # quadind
    append!(instance.l_rows,addedrows)

    rowvec, colvec, coefvec = A.row_pointers, A.columns, A.coefficients
    rows = length(lowerbound)
    cols = size(instance.inner.A)[2]
    push!(rowvec,length(colvec)+1)

    # An = full(sparse(rowvec,colvec,coefvec,rows,cols))
    # @show cols,rows,rowvec,colvec,coefvec
    An = full(SparseMatrixCSC(cols,rows,rowvec,colvec,coefvec)')

    instance.inner.A = vcat(instance.inner.A,An)

    append!(instance.inner.b,lowerbound)
    append!(instance.inner.contype,['R' for i in 1:rows])

    append!(instance.inner.range,upperbound)
    append!(instance.inner.Qcon,[zeros(0,0) for i in 1:rows])
    nothing
end

function modify_ranged_constraints!(instance::MockLinQuadOptimizer, rows, lowerbound, upperbound)
    # quadind
    _rows = instance.l_rows[rows]
    for i in _rows
        instance.inner.b[i] = lowerbound[i]
        instance.inner.range[i] = upperbound[i]
    end
    nothing
end

function LQOI.get_rhs(instance::MockLinQuadOptimizer, row)
    # quadind
    _row = instance.l_rows[row]
    return instance.inner.b[_row]
end
function LQOI.get_quadratic_rhs(instance::MockLinQuadOptimizer, row)
    # quadind
    _row = instance.q_rows[row]
    return instance.inner.b[_row]
end

function LQOI.get_range(instance::MockLinQuadOptimizer, row)
    # quadind
    _row = instance.l_rows[row]
    return instance.inner.b[_row], instance.inner.range[_row]
end

function LQOI.get_linear_constraint(instance::MockLinQuadOptimizer, row)
    # quadind
    _row = instance.l_rows[row]

    vals = instance.inner.A[_row,:]

    outvals = Float64[]
    outinds = Int[]
    for i in eachindex(vals)
        if abs(vals[i]) > 0.0
            push!(outvals,vals[i])
            push!(outinds,i)
        end
    end
    return outinds, outvals
end

function LQOI.get_quadratic_constraint(instance::MockLinQuadOptimizer, row)
    # quadind
    _row = instance.q_rows[row]

    vals = instance.inner.A[_row,:]

    outvals = Float64[]
    outinds = Int[]
    for i in eachindex(vals)
        if abs(vals[i]) > 0.0
            push!(outvals,vals[i])
            push!(outinds,i)
        end
    end

    Q = instance.inner.Qcon[_row]

    #TODO (@joaquim) fix here
    return outinds, outvals, sparse(0.5*Q)
end

function LQOI.change_matrix_coefficient!(instance::MockLinQuadOptimizer, row, col, coef)
    # quadind
    _row = instance.l_rows[row]
    instance.inner.A[_row,col] = coef
end
function LQOI.change_rhs_coefficient!(instance::MockLinQuadOptimizer, row, coef)
    # quadind
    _row = instance.l_rows[row]
    instance.inner.b[_row] = coef
end
function LQOI.change_objective_coefficient!(instance::MockLinQuadOptimizer, col, coef)
    instance.inner.c[col] = coef
end
function LQOI.delete_linear_constraints!(instance::MockLinQuadOptimizer, rowbeg, rowend)
    _rows = instance.l_rows[collect(rowbeg:rowend)]
    survive = setdiff(collect(1:num_cons(instance.inner)), _rows)
    instance.inner.A = instance.inner.A[survive,:]
    instance.inner.b = instance.inner.b[survive]
    instance.inner.contype = instance.inner.contype[survive]
    instance.inner.range = instance.inner.range[survive]
    instance.inner.Qcon = instance.inner.Qcon[survive]
    for i in rowend:-1:rowbeg
        deleteat!(instance.l_rows,i)
    end
end

function LQOI.delete_quadratic_constraints!(instance::MockLinQuadOptimizer, rowbeg, rowend)
    _rows = instance.q_rows[collect(rowbeg:rowend)]
    survive = setdiff(collect(1:num_cons(instance.inner)), _rows)
    instance.inner.A = instance.inner.A[survive,:]
    instance.inner.b = instance.inner.b[survive]
    instance.inner.contype = instance.inner.contype[survive]
    instance.inner.range = instance.inner.range[survive]
    instance.inner.Qcon = instance.inner.Qcon[survive]
    for i in rowend:-1:rowbeg
        deleteat!(instance.q_rows,i)
    end
end

# TODO fix types
function variabletype(::MockLinQuadOptimizer, typeval)
    if typeval == 'B'
        return 'B'
    elseif typeval == 'I'
        return 'I'
    else
        return 'C'
    end
end
function LQOI.change_variable_types!(instance::MockLinQuadOptimizer, colvec, typevec)
    for i in eachindex(colvec)
        instance.inner.vartype[colvec[i]] = variabletype(instance, typevec[i])
    end
end

function sensetype(::MockLinQuadOptimizer, typeval)
    if typeval == 'G'
        return 'G'
    elseif typeval == 'E'
        return 'E'
    elseif typeval == 'L'
        return 'L'
    else
        return 'R'
    end
end
function LQOI.change_linear_constraint_sense!(instance::MockLinQuadOptimizer, rowvec, sensevec)
    _rows = instance.l_rows[rowvec]
    for i in eachindex(_rows)
        instance.inner.vartype[_rows[i]] = sensetype(instance, sensevec[i])
    end
end

function LQOI.add_sos_constraint!(instance::MockLinQuadOptimizer, colvec, valvec, typ)
    sos = LinQuadSOS(typ, valvec, colvec)
    push!(instance.inner.sos,sos)
end
function LQOI.delete_sos!(instance::MockLinQuadOptimizer, idx1, idx2)
    delete_between(instance.inner.sos, idx1, idx2)
end
function LQOI.get_sos_constraint(instance::MockLinQuadOptimizer, idx)
    sos = instance.inner.sos[idx]
    return sos.indices, sos.weights, sos.kind
end

function LQOI.add_quadratic_constraint!(instance::MockLinQuadOptimizer, cols, coefs, rhs, sense, I, J, V)
    @assert length(I) == length(J) == length(V)

    nvars = length(instance.inner.c)
    Q = full(sparse(I,J,V,nvars,nvars))

    a = zeros(nvars)
    for i in eachindex(cols)
        a[cols[i]] = coefs[i]
    end

    instance.inner.A = vcat(instance.inner.A,a')
    push!(instance.inner.b,rhs)
    push!(instance.inner.contype,sense)
    push!(instance.inner.range,Inf)
    push!(instance.inner.Qcon,Q)

    # scalediagonal!(V, I, J, 0.5)
    # scalediagonal!(V, I, J, 2.0)
    push!(instance.q_rows, num_cons(instance.inner))
end

#=
    Objective
=#

function LQOI.set_quadratic_objective!(instance::MockLinQuadOptimizer, I, J, V)
    @assert length(I) == length(J) == length(V)
    n = num_vars(instance.inner)
    Q = full(sparse(I,J,V,n,n))
    # scalediagonal!(V, I, J, 0.5)
    instance.inner.Qobj = Q
    # scalediagonal!(V, I, J, 2.0)
    return nothing
end

function LQOI.set_linear_objective!(instance::MockLinQuadOptimizer, colvec, coefvec)
    nvars = num_vars(instance.inner)
    obj = zeros(Float64, nvars)

    for i in eachindex(colvec)
        obj[colvec[i]] = coefvec[i]
    end

    instance.inner.c = obj
    nothing
end

function LQOI.change_objective_sense!(instance::MockLinQuadOptimizer, symbol)
    if symbol == :min
        instance.inner.sense = :minimize
    else
        instance.inner.sense = :maximize
    end
end

function LQOI.get_linear_objective!(instance::MockLinQuadOptimizer, x)
    for i in eachindex(x)
        x[i] = instance.inner.c[i]
    end
end

function LQOI.get_quadratic_terms_objective(instance::MockLinQuadOptimizer)
    I = Int32[]
    J = Int32[]
    V = Float64[]
    Q = instance.inner.Qobj
    n, m = size(Q)
    for i in 1:n, j in 1:n
        if Q[i,j] != 0
            push!(I,i)
            push!(J,j)
            push!(V,Q[i,j])
        end
    end
    return sparse(Q)
end

function LQOI.get_objectivesense(instance::MockLinQuadOptimizer)
    s = instance.inner.sense
    if s == :maximize
        return MOI.MaxSense
    else
        return MOI.MinSense
    end
end

#=
    Variables
=#

LQOI.get_number_variables(instance::MockLinQuadOptimizer) = length(instance.inner.c)

function LQOI.add_variables!(instance::MockLinQuadOptimizer, int)
    append!(instance.inner.ub,Inf*ones(int))
    append!(instance.inner.lb,-Inf*ones(int))
    append!(instance.inner.c,0.0*ones(int))
    append!(instance.inner.vartype,['C' for i in 1:int])
    instance.inner.A = hcat(instance.inner.A,zeros(length(instance.inner.b),int))
end

function LQOI.delete_variables!(instance::MockLinQuadOptimizer, colbeg, colend)
    instance.inner.Qobj = delete_colrow(instance.inner.Qobj, colbeg, colend)
    instance.inner.Qcon = map(x->delete_colrow(x, colbeg, colend), instance.inner.Qcon)
    if colbeg > 1 && colend < length(instance.inner.lb)
        instance.inner.A = hcat(instance.inner.A[:,1:colbeg-1],instance.inner.A[:,colend+1:end])
        instance.inner.lb = vcat(instance.inner.lb[1:colbeg-1],instance.inner.lb[colend+1:end])
        instance.inner.ub = vcat(instance.inner.ub[1:colbeg-1],instance.inner.ub[colend+1:end])
        instance.inner.c = vcat(instance.inner.c[1:colbeg-1],instance.inner.c[colend+1:end])
        instance.inner.vartype = vcat(instance.inner.vartype[1:colbeg-1],instance.inner.vartype[colend+1:end])
    elseif colbeg > 1
        instance.inner.A = instance.inner.A[:,1:colbeg-1]
        instance.inner.lb = instance.inner.lb[1:colbeg-1]
        instance.inner.ub = instance.inner.ub[1:colbeg-1]
        instance.inner.c = instance.inner.c[1:colbeg-1]
        instance.inner.vartype = instance.inner.vartype[1:colbeg-1]
    elseif colend < length(instance.inner.lb)
        instance.inner.A = instance.inner.A[:,colend+1:end]
        instance.inner.lb = instance.inner.lb[colend+1:end]
        instance.inner.ub = instance.inner.ub[colend+1:end]
        instance.inner.c = instance.inner.c[colend+1:end]
        instance.inner.vartype = instance.inner.vartype[colend+1:end]
    else
        instance.inner.A = zeros(length(instance.inner.b),0)
        instance.inner.lb = Vector{Float64}[]
        instance.inner.ub = Vector{Float64}[]
        instance.inner.c = Vector{Float64}[]
        instance.inner.vartype = Cchar[]
    end
end

function delete_colrow(Q, colbeg, colend)
    n, m = size(Q)
    @assert n == m

    if colbeg > n
        # do nothing
        return Q
    elseif colbeg > 1 && colend < n
        Q = hcat(Q[:,1:colbeg-1], Q[:,colend+1:end])
        Q = vcat(Q[1:colbeg-1,:], Q[colend+1:end,:])
        return Q
    elseif colbeg > 1
        return Q[1:colbeg-1, 1:colbeg-1]
    elseif colend < n
        return Q[colend+1:end, colend+1:end]
    else
        return zeros(0,0)
    end
end

function delete_between(V, colbeg, colend)
    n = length(V)
    if colbeg > n
        # do nothing
        return V
    elseif colbeg > 1 && colend < n
        V = vcat(V[1:colbeg-1], V[colend+1:end])
        return V
    elseif colbeg > 1
        return V[1:colbeg-1]
    elseif colend < n
        return V[colend+1:end]
    else
        return zeros(0)
    end
end

#=
    Solve
=#

LQOI.solve_mip_problem!(instance::MockLinQuadOptimizer) = nothing # fakesolve(instance.inner)

LQOI.solve_quadratic_problem!(instance::MockLinQuadOptimizer) = nothing # fakesolve(instance.inner)

LQOI.solve_linear_problem!(instance::MockLinQuadOptimizer) = nothing # fakesolve(instance.inner)

function LQOI.get_termination_status(instance::MockLinQuadOptimizer)
    return MOI.Success
end

function LQOI.get_primal_status(instance::MockLinQuadOptimizer)
    return MOI.FeasiblePoint
end

function LQOI.get_dual_status(instance::MockLinQuadOptimizer)
    return MOI.FeasiblePoint
end

function LQOI.get_variable_primal_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.variable_primal_solution[i]
    end
    nothing
end

function LQOI.get_linear_primal_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.constraint_primal_solution[i]
    end
    nothing
end

function LQOI.get_quadratic_primal_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.quadratic_primal_solution[i]
    end
    nothing
end

function LQOI.get_variable_dual_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.variable_dual_solution[i]
    end
    nothing
end

function LQOI.get_linear_dual_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.constraint_dual_solution[i]
    end
    nothing
end

function LQOI.get_quadratic_dual_solution!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.quadratic_dual_solution[i]
    end
    nothing
end

LQOI.get_objective_value(instance::MockLinQuadOptimizer) = get_objval(instance.inner)

function LQOI.get_farkas_dual!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.ray_dual_solution[i]
    end
    nothing
end

function hasdualray(instance::MockLinQuadOptimizer)
    return !isempty(instance.inner.ray_dual_solution)
end

function LQOI.get_unbounded_ray!(instance::MockLinQuadOptimizer, place)
    for i in eachindex(place)
        place[i] = instance.inner.ray_primal_solution[i]
    end
    nothing
end

function hasprimalray(instance::MockLinQuadOptimizer)
    return !isempty(instance.inner.ray_primal_solution)
end

MOI.free!(m::MockLinQuadOptimizer) = nothing

#=
LQOI.get_objective_bound(instance::MockLinQuadOptimizer) = get_objval(instance.inner)

function LQOI.get_relative_mip_gap(instance::MockLinQuadOptimizer)
    L = get_objval(instance.inner)
    U = get_objbound(instance.inner)
    return abs(U-L)/U
end

LQOI.get_iteration_count(instance::MockLinQuadOptimizer)  = get_iter_count(instance.inner)

LQOI.get_barrier_iterations(instance::MockLinQuadOptimizer) = get_barrier_iter_count(instance.inner)

LQOI.get_node_count(instance::MockLinQuadOptimizer) = get_node_count(instance.inner)
=#
#=
function LQOI.add_mip_starts!(instance::MockLinQuadOptimizer, colvec::Vector{Int}, valvec::Vector)
    x = zeros(num_vars(instance.inner))
    for (col, val) in zip(colvec, valvec)
        x[col] = val
    end
    loadbasis(instance.inner, x)
end
=#


text:
abstract type InitializationStrategyGamma <: InitializationStrategy end

struct FixedGammaStrategy <: InitializationStrategyGamma
    kernel::SquaredExponentialKernel
end

calculate_gamma(model, strategy::FixedGammaStrategy) = MLKernels.getvalue(strategy.kernel.alpha)

"""
Original publication:
Silverman, Bernard W. Density estimation for statistics and data analysis. Routledge, 2018.
"""
struct RuleOfThumbSilverman <: InitializationStrategyGamma end

function calculate_gamma(model, strategy::RuleOfThumbSilverman)
    return (size(model.data, 2) * (size(model.data, 1) + 2) / 4.0)^(-1.0 / (size(model.data,1) + 4.0))
end

"""
Original publication:
Scott, David W. Multivariate density estimation: theory, practice, and visualization. John Wiley & Sons, 2015.
"""
struct RuleOfThumbScott <: InitializationStrategyGamma end

function calculate_gamma(model, strategy::RuleOfThumbScott)
    return size(model.data, 2)^(-1.0/(size(model.data,1) + 4))
end


text:
@testset "Prints" begin
    random_series = rand(100)
    # Default show methods
    model = LocalLevel(random_series)
    show(model)
    println("")
    fit!(model)
    show(results(model))
    println("")
    # Specific show methods
    model = SARIMA(random_series)
    show(model)
end


text:
# extensions.jl

abstract Extension

immutable PeriodicExtension <: Extension
end

immutable ZeroPadding <: Extension
end

immutable ConstantPadding{T} <: Extension
    constant    ::  T
end

immutable SymmetricExtension{SYM_LEFT,SYM_RIGHT,PAR_LEFT,PAR_RIGHT} <: Extension
end

left_parity{SL,SR,PL,PR}(s::SymmetricExtension{SL,SR,PL,PR})    = PL
right_parity{SL,SR,PL,PR}(s::SymmetricExtension{SL,SR,PL,PR})   = PR
left_symmetry{SL,SR,PL,PR}(s::SymmetricExtension{SL,SR,PL,PR})  = SL
right_symmetry{SL,SR,PL,PR}(s::SymmetricExtension{SL,SR,PL,PR}) = SR

first_index(a) = first(linearindices(a))
last_index(a) = last(linearindices(a))

"""
Return element `k` from the extension sequence induced by `a`.
"""
element(s::Extension, a, k) = _element(s, a, k, first_index(a), last_index(a))

# Note that we do a boundscheck for `k` here. If `k` is a valid index for `a` we
# use it right away, avoiding other (more expensive) index calculations.
_element(s::Extension, a, k, i0, i1) = (i0 <= k <= i1) ? a[k] : element_extension(s, a, k, i0, i1)

# For periodic functions, extension amounts to a modulo operation.
periodic_index(k, i0, i1) = mod(k-i0, i1-i0+1)+i0

element_extension(s::PeriodicExtension, a, k, i0, i1) = a[periodic_index(k, i0, i1)]

# For zero padding of course we return zero (of the same type as elements of a)
element_extension(s::ZeroPadding, a, k, i0, i1) = zero(eltype(a))

# For constant padding we return the stored constant, converted to eltype(a)
element_extension(s::ConstantPadding, a, k, i0, i1) = convert(eltype(a), s.constant)

# The symmetric case is a little more involved.
#
# Methodology: wenever the index `k` is outside the allowed range, we fold it back
# using symmetry and we iterate. With each iteration, it goes closer to the allowed range.
# An alternative for large indices could be based on the periodicity of the symmetric
# extension.

function _element(s::SymmetricExtension, a, k, i0, i1)
    if k > i1
        # We are to the right of the range i0:i1
        element_right(s, a, k, i0, i1)
    elseif k < i0
        # We are to the left of the range i0:i1
        element_left(s, a, k, i0, i1)
    else
        # k is a valid index for a
        a[k]
    end
end

# Right whole point symmetry, even or odd
element_right{SYM_LEFT,PAR_LEFT}(s::SymmetricExtension{SYM_LEFT,:wp,PAR_LEFT,:even}, a, k, i0, i1) = _element(s, a, 2*i1 - k + 1, i0, i1)
element_right{SYM_LEFT,PAR_LEFT}(s::SymmetricExtension{SYM_LEFT,:wp,PAR_LEFT,:odd}, a, k, i0, i1) = -_element(s, a, 2*i1 - k + 1, i0, i1)

# Right half point symmetry, even or odd
element_right{SYM_LEFT,PAR_LEFT}(s::SymmetricExtension{SYM_LEFT,:hp,PAR_LEFT,:even}, a, k, i0, i1) = _element(s, a, 2*i1 - k, i0, i1)
element_right{SYM_LEFT,PAR_LEFT}(s::SymmetricExtension{SYM_LEFT,:hp,PAR_LEFT,:odd}, a, k, i0, i1) = -_element(s, a, 2*i1 - k, i0, i1)

# Left whole point symmetry, even or odd
element_left{SYM_RIGHT}(s::SymmetricExtension{:wp,SYM_RIGHT,:even}, a, k, i0, i1) = _element(s, a, 2*i0 - k - 1, i0, i1)
element_left{SYM_RIGHT}(s::SymmetricExtension{:wp,SYM_RIGHT,:odd}, a, k, i0, i1) = -_element(s, a, 2*i0 - k - 1, i0, i1)

# Left half point symmetry, even or odd
element_left{SYM_RIGHT}(s::SymmetricExtension{:hp,SYM_RIGHT,:even}, a, k, i0, i1) = _element(s, a, 2*i0 - k, i0, i1)
element_left{SYM_RIGHT}(s::SymmetricExtension{:hp,SYM_RIGHT,:odd}, a, k, i0, i1) = -_element(s, a, 2*i0 - k, i0, i1)


text:
# TO RUN:
# include("main.jl")
#
# AUTHOR
# Dr. Umberto Saetti, Assistant Professor, Department of Aerospace Engineering,
# Auburn university
# saetti@auburn.edu
#
# LAST UPDATED
# 1/3/2022
#
# DESCRIPTION
# Nonlinear simulation model of the rigid-body dynamics of a conventional
# helicopter configuration (i.e., main rotor + tail rotor).
# The model uses static aerodynamic models for fuselage, main rotor, tail
# rotor, and empennage. The trim, linearization routines, and overall
# architecture of the simulation follows the teachings of Dr. Joe Horn at Penn
# State. A more detailed description of the model, including references can be
# found in Ref. 1. If using the code, please cite Ref. 1.
#
# REFERENCES
# 1) Saetti, U., and Horn, J. F., "Flight Simulation and Control using the Julia
#    Language", AIAA Scitech Forum, San Diego, CA, Jan 3-7, 2022.
#    DOI: https://arc.aiaa.org/doi/10.2514/6.2022-2354.
#
# STATES             | NINDICES    | DESCRIPTION
# ______________________________________________________________________________
# u v w              |  1  2  3    | body velocities [ft/s]
# p q r              |  4  5  6    | angula rates [rad]
# phi theta psi      |  7  8  9    | Euler angles [rad]
# x  y  z            | 10 11 12    | position [ft]
#
# CONTROL INPUTS     | INDICED     | DESCRIPTION
# ______________________________________________________________________________
# delta_lat          | 1           | lateral stick [%]
# delta_lon          | 2           | longitudinal stick [%]
# delta_col          | 3           | collective stick [%]
# delta_ped          | 4           | pedals [%]
#
# ------------------------------------------------------------------------------

# include packages
#using DifferentialEquations
using LinearAlgebra
using SparseArrays
using ControlSystems
using MAT
using Interpolations
using ApproxFun
using Printf
#ENV["MPLBACKEND"]="tkagg"
using Plots
using JLD
#using PyPlot
#using GR

@printf("\n                 J-SimpleHel           \n")
@printf("\n          Author: Dr. Umberto Saetti\n")
@printf("\n-----------------------------------------------\n")

# include functions
include("atan2.jl")
include("linearize.jl")
include("trimmer.jl")
include("SimpleHel.jl")
include("control_mixing.jl")
include("SimpleHel_DI.jl")
include("rk4.jl")
include("finp.jl")
include("simulate.jl")
include("simulate_DI.jl")
include("wrapper.jl")
# load H-60 constants
include("H60_constants.jl")

# ------------------------ TRIM FLIGHT DYNAMICS MODEL --------------------------

# trim speeds [ft/s]
VXTRIM=80*1.688
VYTRIM=0.0
VZTRIM=0.0
# trim turn rate [rad/s]
PSIDTRIM=0.0
include("trimInit.jl")
# trim aircraft
x0, u0 = trimmer(SimpleHel!,x0,u0,xdot_targ)
# linearize aircraft
A, B, C, D = linearize(SimpleHel!,x0,u0)
# partition in longitudinal and lateral flight dynamics
Alon=A[[1,3,5,8],[1,3,5,8]]
Alat=A[[2,4,6,7],[2,4,6,7]]

# ----------------------------- SPECTRAL ANALYSIS ------------------------------

# eigenvalues
eigs=eigvals(A)
eigsLon=eigvals(Alon)
eigsLat=eigvals(Alat)
# store eigenvalues
#save("AB_J-SimpleHel.jld","A",A,"B",B)
#save("eigs_J-SimpleHel.jld","eigs",eigs)

# plot eigenvalues
gr(size=(800,600))
plt_eigs=plot(real(eigs),imag(eigs),seriestype=:scatter,
      label="Full-Order", legendfontsize=12, legend=:topleft,
      marker = (:circle, :royalblue, 6), markerstrokecolor = :royalblue)
plot!(real(eigsLat),imag(eigsLat),seriestype=:scatter,
      label="Reduced-Order (Lat)", legendfontsize=12,
      marker = (:utriangle, :brown3, 6), markerstrokecolor = :brown3,
      reuse = false, color=:purple)
plot!(real(eigsLon),imag(eigsLon),seriestype=:scatter,
      label="Reduced-Order (Lon)", legendfontsize=12,
      marker = (:star5, :forestgreen, 6), markerstrokecolor = :forestgreen)
display(plt_eigs)
xaxis!("Real",xguidefontsize=14)
yaxis!("Imag",yguidefontsize=14)

# ---------------------------- FREQUENCY RESPONSES -----------------------------

# linear system
sys=ss(A,B,Matrix{Float64}(I, NSTATES, NSTATES),zeros(NSTATES,NCTRLS))
# transfer function (p/dlat)
pdlat=tf(sys[4,1])
# freq response
mag, phase, w = bode(pdlat, 10.0.^(range(-1,stop=2,length=1000)))

# plot pdlat
gr(size=(1000,1000))
plt_mag=plot(vec(w), 20*log10.(vec(mag)/perc2in[1,1]), label="Full order",
             reuse = false, line=(4, :black, :solid))
xaxis!(:log10)
yaxis!("Magnitude [dB]")
if VXTRIM==80*1.688
      title!("p/dlat, KTAS = 80 kts")
elseif VXTRIM<=1*1.688
      title!("p/dlat, KTAS = 0 kts")
end
xlims!((0.1,100))
plt_phase=plot(vec(w), wrapper(vec(phase)).-360, label="", line=(4, :black, :solid))
yaxis!("Phase [deg]")
xaxis!(:log10)
xlims!((0.1,100))
plt_pdlat=plot(plt_mag, plt_phase, layout=(2,1), margin=10Plots.mm)
display(plt_pdlat)
#if VXTRIM==80*1.688
#      savefig(".\\Plots\\pdlat_80kt.png")
#elseif VXTRIM<=1*1.688
#      savefig(".\\Plots\\pdlat_Hover.png")
#end

# transfer function (q/dlon)
pdlat=tf(sys[5,2])
# freq response
mag, phase, w = bode(pdlat, 10.0.^(range(-1,stop=2,length=1000)))

# plot qdlon
gr(size=(1000,1000))
plt_mag=plot(vec(w), 20*log10.(vec(mag)/perc2in[2,2]), label="Full order",
             reuse = false, line=(4, :black, :solid))
xaxis!(:log10)
yaxis!("Magnitude [dB]")
if VXTRIM==80*1.688
      title!("q/dlon, KTAS = 80 kts")
elseif VXTRIM<=1*1.688
      title!("q/dlon, KTAS = 0 kts")
end
xlims!((0.1,100))
plt_phase=plot(vec(w), wrapper(vec(phase)).-360, label="", line=(4, :black, :solid))
yaxis!("Phase [deg]")
xaxis!(:log10)
xlims!((0.1,100))
plt_qdlon=plot(plt_mag, plt_phase, layout=(2,1), margin=10Plots.mm)
display(plt_qdlon)
#if VXTRIM==80*1.688
#      savefig(".\\Plots\\qdlon_80kt.png")
#elseif VXTRIM<=1*1.688
#      savefig(".\\Plots\\qdlon_Hover.png")
#end

# ------------------------------ TIME SIMULATION -------------------------------

# time step [s]
dt=0.01
# length of simulation [s]
Tsim=10
# run open-loop simulation
state_OL, time_OL = simulate(SimpleHel!,finp,Tsim,dt,x0,u0)
# number of dynamic inverse states
NDISTATES=9
# initial state vector of closed-loop simulation
x0_DI=[x0; 0.0; 0.0; 0.0; 0.0; 0.0; 0.0; u0[1]; u0[2]; u0[4]]
# run closed-loop simulation
state_CL, time_CL = simulate_DI(SimpleHel_DI!,finp,Tsim,dt,x0,u0)
# store open- and closed-loop simulations
#save("resp_SimpleHel.jld","time_OL",time_OL,"time_CL",time_CL,
#    "state_OL",state_OL,"state_CL",state_CL)

# plot attitude
gr(size=(800,600))
i=7
p1=plot(time_OL,state_OL[i,:]*180/pi,label="Open Loop",legend=:topright, legendfontsize=12,line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:]*180/pi,label="Closed Loop",legend=:topright, legendfontsize=12,line=(3, :brown3, :dash))
yaxis!("φ [deg]",yguidefontsize=14)
xlims!((0,10))
i=8
p2=plot(time_OL,state_OL[i,:]*180/pi,label="",line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:]*180/pi,label="",line=(3, :brown3, :dash))
yaxis!("θ [deg]",yguidefontsize=14)
xlims!((0,10))
i=9
p3=plot(time_OL,state_OL[i,:]*180/pi,label="",line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:]*180/pi,label="",line=(3, :brown3, :dash))
xaxis!("Time [s]",xguidefontsize=14)
yaxis!("ψ [deg]",yguidefontsize=14)
xlims!((0,10))
p4=plot(p1, p2, p3, layout=(3,1), leftmargin=3Plots.mm)
display(p4)
#savefig(".\\Plots\\att_J-SimpleHel.svg")

# plot angular rates
gr(size=(800,600))
i=4
p5=plot(time_OL,state_OL[i,:],label="Open Loop",legend=:bottomleft, legendfontsize=12,line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:],label="Closed Loop",legend=:bottomleft, legendfontsize=12,line=(3, :brown3, :dash))
yaxis!("p [rad/s]",yguidefontsize=14)
xlims!((0,10))
i=5
p6=plot(time_OL,state_OL[i,:],label="",line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:],label="",line=(3, :brown3, :dash))
yaxis!("q [rad/s]",yguidefontsize=14)
xlims!((0,10))
i=6
p7=plot(time_OL,state_OL[i,:],label="",line=(3, :royalblue, :solid))
plot!(time_OL,state_CL[i,:],label="",line=(3, :brown3, :dash))
xaxis!("Time [s]",xguidefontsize=14)
yaxis!("r [rad/s]",yguidefontsize=14)
xlims!((0,10))
p8=plot(p5, p6, p7, layout=(3,1), leftmargin=3Plots.mm)
display(p8)
#savefig(".\\Plots\\ang_J-SimpleHel.svg")

max_stars_repo_name: umbertosaetti/J-SimpleHel

max_issues_repo_name: umbertosaetti/J-SimpleHel

max_forks_repo_name: umbertosaetti/J-SimpleHel


text:
export Spindle

"""
`Spindle()` returns the Moser spindle graph. This is a seven-vertex
unit distance graph with chromatic number equal to 4.
"""
function Spindle()
    G = IntGraph(7)
    edges = [
        1 2
        1 3
        2 3
        2 4
        3 4
        1 5
        1 6
        5 6
        5 7
        6 7
        4 7 ]
    add_edges!(G,edges)

    d = Dict{Int,Vector}()
    a = sqrt(3)/2

    pts = [ 0 1/2 -1/2 0 ; 0 a a 2a ]

    theta = acos(5/6)/2
    R = [ cos(theta) -sin(theta); sin(theta) cos(theta) ]

    p1 = R*pts
    for k=1:4
        d[k] = p1[:,k]
    end

    p2 = R'*pts
    for k=5:7
        d[k] = p2[:,k-3]
    end
    embed(G,d)
    SimpleGraphs.name(G,"Moser Spindle")
    return G
end


text:
function ispoweroftwo(n)
  if n > 1
    while n & 1 == 0
      n >>= 1
    end
  end
  n == 1
end

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
"""
    ft_center_diff(s [, dims])

Return the shifts required to center dimensions `dims` at the respective
Fourier centers.
This function is internally used by [`ShiftedArrays.fftshift`](@ref) and
[`ShiftedArrays.ifftshift`](@ref).

# Examples

```jldoctest
julia> ShiftedArrays.ft_center_diff((4, 5, 6), (1, 2)) # Fourier center is at (2, 3, 0)
(2, 2, 0)

julia> ShiftedArrays.ft_center_diff((4, 5, 6), (1, 2, 3)) # Fourier center is at (2, 3, 4)
(2, 2, 3)
```
"""
function ft_center_diff(s::NTuple{N, T}, dims=ntuple(identity, Val(N))) where {N, T}
    return ntuple(i -> i ∈ dims ?  s[i] ÷ 2 : 0, N)
end

"""
    fftshift(x [, dims])

Lazy version of `AbstractFFTs.fftshift(x, dims)`. Return a `CircShiftedArray`
where each given dimension is shifted by `N÷2`, where `N` is the size of
that dimension.

# Examples

```jldoctest
julia> ShiftedArrays.fftshift([1 0 0 0])
1×4 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 0  0  1  0

julia> ShiftedArrays.fftshift([1 0 0; 0 0 0; 0 0 0])
3×3 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 0  0  0
 0  1  0
 0  0  0

julia> ShiftedArrays.fftshift([1 0 0; 0 0 0; 0 0 0], (1,))
3×3 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 0  0  0
 1  0  0
 0  0  0
```
"""
function fftshift(x::AbstractArray{T, N}, dims=ntuple(identity, Val(N))) where {T, N}
    ShiftedArrays.circshift(x, ft_center_diff(size(x), dims))
end

"""
    ifftshift(x [, dims])

Lazy version of `AbstractFFTs.ifftshift(x, dims)`. Return a `CircShiftedArray`
where each given dimension is shifted by `-N÷2`, where `N` is the size of
that dimension.

# Examples

```jldoctest
julia> ShiftedArrays.ifftshift([0 0 1 0])
1×4 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 1  0  0  0

julia> ShiftedArrays.ifftshift([0 0 0; 0 1 0; 0 0 0])
3×3 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 1  0  0
 0  0  0
 0  0  0

julia> ShiftedArrays.ifftshift([0 1 0; 0 0 0; 0 0 0], (2,))
3×3 CircShiftedArray{Int64, 2, Matrix{Int64}}:
 1  0  0
 0  0  0
 0  0  0
```
"""
function ifftshift(x::AbstractArray{T, N}, dims=ntuple(identity, Val(N))) where {T, N}
    ShiftedArrays.circshift(x, map(-, ft_center_diff(size(x), dims)))
end


text:
function reset!(ai::AbstractInitial) end

function update!(ai::AbstractInitial, state, p, α) end

function reset!(alg::LineSearch)
    reset!(alg.initial)
end

function update!(alg, state, p, α)
    update!(alg.initial, state, p, α)
    return nothing
end

function step!(x, alg, f::Function, state)
    p = direction(alg.direction, state)
    α₀ = alg.initial(state, p)
    α = alg.method(f, state, p, α₀)
    update!(alg, state, p, α)
    x .+= α*p
end


text:
using AtomsBase
# Key functionality to integrate DFTK and AtomsBase

function construct_system(lattice::AbstractMatrix, atoms::Vector, positions::Vector, magnetic_moments::Vector)
    @assert length(atoms) == length(positions)
    @assert isempty(magnetic_moments) || length(magnetic_moments) == length(atoms)
    atomsbase_atoms = map(enumerate(atoms)) do (i, element)
        kwargs = Dict{Symbol, Any}(:potential => element, )
        if !isempty(magnetic_moments)
            kwargs[:magnetic_moment] = normalize_magnetic_moment(magnetic_moments[i])
        end
        if element isa ElementPsp
            kwargs[:pseudopotential] = element.psp.identifier
        end

        position = lattice * positions[i] * u"bohr"
        if atomic_symbol(element) == :X  # dummy element ... should solve this upstream
            Atom(:X, position; atomic_symbol=:X, atomic_number=0, atomic_mass=0u"u", kwargs...)
        else
            Atom(atomic_symbol(element), position; kwargs...)
        end
    end
    periodic_system(atomsbase_atoms, collect(eachcol(lattice)) * u"bohr")
end

function parse_system(system::AbstractSystem{D}) where {D}
    if !all(periodicity(system))
        error("DFTK only supports calculations with periodic boundary conditions.")
    end

    # Parse abstract system and return data required to construct model
    mtx = austrip.(reduce(hcat, bounding_box(system)))
    T = eltype(mtx)
    lattice = zeros(T, 3, 3)
    lattice[1:D, 1:D] .= mtx

    # Cache for instantiated pseudopotentials. This is done to ensure that identical
    # atoms are indistinguishable in memory, which is used in the Model constructor
    # to deduce the atom_groups.
    cached_pspelements = Dict{String, ElementPsp}()
    atoms = map(system) do atom
        if hasproperty(atom, :potential)
            atom.potential
        elseif hasproperty(atom, :pseudopotential)
            get!(cached_pspelements, atom.pseudopotential) do
                ElementPsp(atomic_symbol(atom); psp=load_psp(atom.pseudopotential))
            end
        else
            ElementCoulomb(atomic_symbol(atom))
        end
    end

    positions = map(system) do atom
        coordinate = zeros(T, 3)
        coordinate[1:D] = lattice[1:D, 1:D] \ T.(austrip.(position(atom)))
        Vec3{T}(coordinate)
    end

    magnetic_moments = map(system) do atom
        hasproperty(atom, :magnetic_moment) || return nothing
        getproperty(atom, :magnetic_moment)
    end
    if all(m -> isnothing(m) || iszero(m) || isempty(m), magnetic_moments)
        empty!(magnetic_moments)
    else
        magnetic_moments = normalize_magnetic_moment.(magnetic_moments)
    end

    # TODO Use system to determine n_electrons

    (; lattice, atoms, positions, magnetic_moments)
end

function _call_with_system(f, system::AbstractSystem, args...; kwargs...)
    @assert !(:magnetic_moments in keys(kwargs))
    parsed = parse_system(system)
    f(parsed.lattice, parsed.atoms, parsed.positions, args...;
      parsed.magnetic_moments, kwargs...)
end

"""
    Model(system::AbstractSystem; kwargs...)

AtomsBase-compatible Model constructor. Sets structural information (`atoms`, `positions`,
`lattice`, `n_electrons` etc.) from the passed `system`.
"""
Model(system::AbstractSystem; kwargs...) = _call_with_system(Model, system; kwargs...)

# Generate equivalent functions for AtomsBase
for fun in (:model_atomic, :model_DFT, :model_LDA, :model_PBE, :model_SCAN)
    @eval function $fun(system::AbstractSystem, args...; kwargs...)
        _call_with_system($fun, system, args...; kwargs...)
    end
end

# Extra methods to AtomsBase functions for DFTK data structures
"""
    atomic_system(model::DFTK.Model, magnetic_moments=[])

Construct an AtomsBase atomic system from a DFTK model and associated magnetic moments.
"""
function AtomsBase.atomic_system(model::Model, magnetic_moments=[])
    construct_system(model.lattice, model.atoms, model.positions, magnetic_moments)
end

"""
    periodic_system(model::DFTK.Model, magnetic_moments=[])

Construct an AtomsBase atomic system from a DFTK model and associated magnetic moments.
"""
function AtomsBase.periodic_system(model::Model, magnetic_moments=[])
    atomic_system(model, magnetic_moments)
end

AtomsBase.chemical_formula(model::Model) = chemical_formula(atomic_symbol.(model.atoms))


text:
@testset "PolyVar order" begin
    @polyvar x
    z = x
    @polyvar x
    @test z != x
end
@testset "Same variables grlex" begin
    @test DynamicPolynomials.samevars_grlex([1, 0, 1], [1, 1, 0]) == -1
    @test DynamicPolynomials.samevars_grlex([1, 1, 0], [1, 0, 1]) == 1
end


text:
using SQLite
using NeXLCore
using IntervalSets

struct DBFitSpectrum
    database::SQLite.DB
    campaign::Int
    spectrum::DBSpectrum
end

Base.show(io::IO, dbf::DBFitSpectrum) = print(io, "Fit[$(dbf.spectrum.name)]")

NeXLUncertainties.asa(::Type{Spectrum}, dbfs::DBFitSpectrum) = NeXLUncertainties.asa(Spectrum, dbfs.spectrum)

function NeXLUncertainties.asa(::Type{DataFrame}, refs::Vector{DBFitSpectrum})
    fm, sm, cm = Union{Float64, Missing}, Union{String, Missing}, Union{Material,Missing}
    spname, e0, lt, pc, comp = sm[], fm[],  fm[], fm[], cm[]
    for ref in refs
        spec = asa(Spectrum, ref)
        push!(spname, get(spec, :Name, missing))
        push!(e0, get(spec, :BeamEnergy, missing))
        push!(lt, get(spec, :LiveTime, missing))
        push!(pc, get(spec, :ProbeCurrent, missing))
        push!(comp, get(spec, :Composition, missing))
    end
    return DataFrame(Name=spname, BeamEnergy=e0, LiveTime=lt, ProbeCurrent=pc, Composition=comp)
end

function Base.write(db::SQLite.DB, ::Type{DBFitSpectrum}, campaign::Int, spectrum::Int)::Int
    stmt1 = SQLite.Stmt(db, "INSERT INTO FITSPECTRUM(CAMPAIGN, SPECTRUM ) VALUES ( ?, ? );")
    q = DBInterface.execute(stmt1, (campaign, spectrum))
    return DBInterface.lastrowid(q)
end


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

mutable struct IoK8sApiCoreV1GCEPersistentDiskVolumeSource <: SwaggerModel
    fsType::Any # spec type: Union{ Nothing, String } # spec name: fsType
    partition::Any # spec type: Union{ Nothing, Int32 } # spec name: partition
    pdName::Any # spec type: Union{ Nothing, String } # spec name: pdName
    readOnly::Any # spec type: Union{ Nothing, Bool } # spec name: readOnly

    function IoK8sApiCoreV1GCEPersistentDiskVolumeSource(;fsType=nothing, partition=nothing, pdName=nothing, readOnly=nothing)
        o = new()
        validate_property(IoK8sApiCoreV1GCEPersistentDiskVolumeSource, Symbol("fsType"), fsType)
        setfield!(o, Symbol("fsType"), fsType)
        validate_property(IoK8sApiCoreV1GCEPersistentDiskVolumeSource, Symbol("partition"), partition)
        setfield!(o, Symbol("partition"), partition)
        validate_property(IoK8sApiCoreV1GCEPersistentDiskVolumeSource, Symbol("pdName"), pdName)
        setfield!(o, Symbol("pdName"), pdName)
        validate_property(IoK8sApiCoreV1GCEPersistentDiskVolumeSource, Symbol("readOnly"), readOnly)
        setfield!(o, Symbol("readOnly"), readOnly)
        o
    end
end # type IoK8sApiCoreV1GCEPersistentDiskVolumeSource

const _property_map_IoK8sApiCoreV1GCEPersistentDiskVolumeSource = Dict{Symbol,Symbol}(Symbol("fsType")=>Symbol("fsType"), Symbol("partition")=>Symbol("partition"), Symbol("pdName")=>Symbol("pdName"), Symbol("readOnly")=>Symbol("readOnly"))
const _property_types_IoK8sApiCoreV1GCEPersistentDiskVolumeSource = Dict{Symbol,String}(Symbol("fsType")=>"String", Symbol("partition")=>"Int32", Symbol("pdName")=>"String", Symbol("readOnly")=>"Bool")
Base.propertynames(::Type{ IoK8sApiCoreV1GCEPersistentDiskVolumeSource }) = collect(keys(_property_map_IoK8sApiCoreV1GCEPersistentDiskVolumeSource))
Swagger.property_type(::Type{ IoK8sApiCoreV1GCEPersistentDiskVolumeSource }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiCoreV1GCEPersistentDiskVolumeSource[name]))}
Swagger.field_name(::Type{ IoK8sApiCoreV1GCEPersistentDiskVolumeSource }, property_name::Symbol) =  _property_map_IoK8sApiCoreV1GCEPersistentDiskVolumeSource[property_name]

function check_required(o::IoK8sApiCoreV1GCEPersistentDiskVolumeSource)
    (getproperty(o, Symbol("pdName")) === nothing) && (return false)
    true
end

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


text:
"""
        mairca(decisionMat, weights, fns)

Apply MAIRCA (Multi Attributive Ideal-Real Comparative Analysis) for a given matrix and weights.

# Arguments:
 - `decisionMat::DataFrame`: n × m matrix of objective values for n alternatives and m criteria
 - `weights::Array{Float64, 1}`: m-vector of weights that sum up to 1.0. If the sum of weights is not 1.0, it is automatically normalized.
 - `fns::Array{Function, 1}`: m-vector of functions to be applied on the columns.

# Description
mairca() applies the MAIRCA method to rank n alternatives subject to m criteria which are supposed to be
either maximized or minimized.

# Output
- `::MAIRCAResult`: MAIRCAResult object that holds multiple outputs including scores, rankings, and best index.

# Examples
```julia-repl
julia> decmat = [6.952 8.000 6.649 7.268 8.000 7.652 6.316;
       7.319 7.319 6.604 7.319 8.000 7.652 5.313;
       7.000 7.319 7.652 6.952 7.652 6.952 4.642;
       7.319 6.952 6.649 7.319 7.652 6.649 5.000]
4×7 Array{Float64,2}:
 6.952  8.0    6.649  7.268  8.0    7.652  6.316
 7.319  7.319  6.604  7.319  8.0    7.652  5.313
 7.0    7.319  7.652  6.952  7.652  6.952  4.642
 7.319  6.952  6.649  7.319  7.652  6.649  5.0

 julia> df = makeDecisionMatrix(decmat)
 4×7 DataFrame
  Row │ Crt1     Crt2     Crt3     Crt4     Crt5     Crt6     Crt7
      │ Float64  Float64  Float64  Float64  Float64  Float64  Float64
 ─────┼───────────────────────────────────────────────────────────────
    1 │   6.952    8.0      6.649    7.268    8.0      7.652    6.316
    2 │   7.319    7.319    6.604    7.319    8.0      7.652    5.313
    3 │   7.0      7.319    7.652    6.952    7.652    6.952    4.642
    4 │   7.319    6.952    6.649    7.319    7.652    6.649    5.0

julia> weights = [0.172, 0.165, 0.159, 0.129, 0.112, 0.122, 0.140];

julia> fns = [maximum, maximum, maximum, maximum, maximum, maximum, minimum];

julia> result = mairca(df, weights, fns);

julia> result.scores
4-element Array{Float64,1}:
 0.12064543054088471
 0.08066456363291889
 0.14586265389012484
 0.14542366685864686

julia> result.bestIndex
2
```

# References

Pamučar, D., Lukovac, V., Božanić, D., & Komazec, N. (2018). Multi-criteria FUCOM-MAIRCA model for the evaluation of level crossings: case study in the Republic of Serbia. Operational Research in Engineering Sciences: Theory and Applications, 1(1), 108-129.

Ulutaş A.(2019),Swara Ve Mairca Yöntemleri İle Catering Firması Seçimi,BMIJ, (2019), 7(4): 1467-1479 http://dx.doi.org/10.15295/bmij.v7i4.1166
"""
function mairca(decisionMat::DataFrame, weights::Array{Float64,1}, fns::Array{Function,1})::MAIRCAResult

    row, col = size(decisionMat)

    w = unitize(weights)

    T = zeros(Float64, row, col)

    for i in 1:col
        T[:, i] .= w[i] * (1 / row)
    end

    colMax = colmaxs(decisionMat)
    colMin = colmins(decisionMat)

    A = similar(decisionMat)

    for i in 1:row
        for j in 1:col
            if fns[j] == maximum
                @inbounds A[i, j] = T[i,j] * ((decisionMat[i, j] - colMin[j]) / (colMax[j] - colMin[j]))
            elseif fns[j] == minimum
                @inbounds A[i, j] = T[i,j] * ((decisionMat[i, j] - colMax[j]) / (colMin[j] - colMax[j]))
            end
        end
    end

    S = T .- A

    scores = zeros(Float64, row)
    for i in 1:row
        scores[i] = sum(S[i, :])
    end

    rankings = sortperm(scores)

    bestIndex = rankings |> first

    result = MAIRCAResult(
        decisionMat,
        w,
        scores,
        rankings,
        bestIndex
    )

    return result
end

max_stars_repo_name: bahadirfyildirim/JMcDM

max_issues_repo_name: bahadirfyildirim/JMcDM

max_forks_repo_name: bahadirfyildirim/JMcDM


text:
# open a fastq stream, mode should be either r or w
function fastq_open(filename::AbstractString, mode::AbstractString="r")
	return opengene_open(filename, mode)
end

# read the stream and return a fastq record
# \n will be stripped automatically
function fastq_read(stream)
	try
		if eof(stream)
			return false
		end
		name = rstrip(readline(stream),'\n')
		sequence = Sequence(rstrip(readline(stream),'\n'))
		strand = rstrip(readline(stream),'\n')
		quality = Quality(rstrip(readline(stream),'\n'))
		return FastqRead(name, sequence, strand, quality)
	catch e
		println(e)
		return false
	end
end

# write a fastq read to a stream
function fastq_write(stream, fq::FastqRead)
	try
		strs = [fq.name, fq.sequence.seq, fq.strand, fq.quality.qual]
		write(stream, join(strs, "\n")*"\n")
	catch e
		println(e)
		return false
	end
end

# handle reading and writing for a pair<read1 + read2> of fastq files
type FastqPairStream
	read1stream
	read2stream
	iswrite::Bool
end

function fastq_flush_pair(s::FastqPairStream)
	if s.iswrite
		flush(s.read1stream)
		flush(s.read2stream)
	end
end

function fastq_close_pair(s::FastqPairStream)
	if s.iswrite
		close(s.read1stream)
		close(s.read2stream)
	end
end

function fastq_open_pair(filename1::AbstractString, filename2::AbstractString, mode::AbstractString="r")
	return FastqPairStream(fastq_open(filename1, mode), fastq_open(filename2, mode), contains(mode, "w"))
end

function fastq_read_pair(s::FastqPairStream)
	read1 = fastq_read(s.read1stream)
	read2 = fastq_read(s.read2stream)
	if read1 == false || read2 == false
		return false
	end
	return FastqPair(read1, read2)
end

function fastq_write_pair(s::FastqPairStream, pair::FastqPair)
	fastq_write(s.read1stream, pair.read1)
	fastq_write(s.read2stream, pair.read2)
end

Base.eof(stream::FastqPairStream) =  eof(stream.read1stream) || eof(stream.read2stream)
Base.flush(stream::FastqPairStream) = fastq_flush_pair(stream)
Base.close(stream::FastqPairStream) = fastq_close_pair(stream)


text:
module FindProc

using JuMP
using Cbc

include("Format.jl")
using .Format

export findproc

B = []

# Determines procedure with minimum piece size alpha -- Work in progress
# Author: Jason Liu
function findproc(m::Int64, s::Int64, alpha::Rational{Int64}; output::Int64 = 2)
    # Prepare variables for solving
    c = numerator(alpha)
    d = denominator(alpha)
    b = lcm(s, d)

    a = Int64(b * c / d)
    global B = [x for x = a:(b-a)]

    M = []
    S = []
    procedures = []

    # Display header
    alphaF = formatFrac(alpha)
    if output > 0
        printHeader(center("FIND PROCEDURE"))
        printfT(
            "Overview",
            "Find procedures for dividing $m muffins among $s students where $alphaF is the smallest piece size",
            "",
            "Let the common denominator be $b",
        )
    end

    # Exit early if not solvable via Cbc
    if m <= 0 || s <= 0 || alpha <= 0 || alpha > 1
        if output > 0
            printf("Find Procedure does not apply", line = true)
            printEnd()
        end
        return [b, Nothing]
    elseif m % s == 0 && toFrac(alpha) == 1
        b = 1
        give = repeat([1], Int64(m / s))
        if output > 0
            printfT(
                "Procedure",
                "Divide $m muffins {1}",
                "Give $s students {$(join(give, ","))}",
            )
            printEnd()
        end
        return [b, [[[[m, [1]]], [[s, give]]]]]
    end

    try
        # Determine vector sets
        L = b - 2a + 1
        V = Int64(ceil(2m / s))
        T = Int64(b * m / s)
        M = vectorize(unionF(f(L, b, 2), f(L, b, 3)))
        S = vectorize(unionF(f(L, T, V), f(L, T, V - 1)))

        # Initialize model
        model = Model(Cbc.Optimizer)
        set_optimizer_attribute(model, "logLevel", 0)
        set_optimizer_attribute(model, "allowableGap", 1e-5)
        set_optimizer_attribute(model, "cuts", "off")

        @variable(model, 0 <= x[i = 1:length(M)], Int)
        @variable(model, 0 <= y[i = 1:length(S)], Int)

        for k = 1:L
            @constraint(model, sum(x .* getindex.(M, k)) == sum(y .* getindex.(S, k)))
        end

        @constraint(model, sum(x) == m)
        @constraint(model, sum(y) == s)

        # Solve for a natural number solutions
        for i = 1:min(length(M), 3)
            @objective(model, Min, x[i])
            optimize!(model)
            append!(procedures, [[round.(value.(x)), round.(value.(y))]])

            @objective(model, Max, x[i])
            optimize!(model)
            append!(procedures, [[round.(value.(x)), round.(value.(y))]])
        end

        procedures = unique(procedures)
    catch
        Nothing
    end

    # Exit if no procedures
    if length(procedures) == 0
        if output > 0
            printf("No procedures for muffins($m, $s, $alphaF)", line = true)
            printEnd()
        end
        return [b, Nothing]
    end

    # Output each procedure
    M = [unpack(M[i]) for i = 1:length(M)]
    S = [unpack(S[i]) for i = 1:length(S)]
    formatProc = []
    for k = 1:length(procedures)
        x = procedures[k][1]
        y = procedures[k][2]

        divide = [
            "Divide $(Int64(x[i])) muffins {$(join(M[i], ","))}"
            for i = 1:length(M) if x[i] > 0
        ]
        give = [
            "Give $(Int64(y[i])) students {$(join(S[i], ","))}"
            for i = 1:length(S) if y[i] > 0
        ]
        append!(
            formatProc,
            [[
                [[Int64(x[i]), M[i]] for i = 1:length(M)],
                [[Int64(y[i]), S[i]] for i = 1:length(S)],
            ]],
        )

        if output > 1
            printfT("Procedure $k", divide..., give...)
        elseif output > 0
            printfT("Procedure", divide..., give...)
            printEnd()
            output = 0
        end
    end
    output > 0 && printEnd()

    [b, formatProc]
end

# Helper function for findproc -- "vectorizes" arrays in set S based on B
function vectorize(S)
    S == Nothing ? Nothing :
    sort([[count(i -> (i == x), s) for x in B] for s in S], rev = true)
end

# Helper function for findproc -- "unpacks" vectors into strings based on B
function unpack(V)
    vcat(map(i -> repeat([B[i]], V[i]), 1:length(V))...)
end

# Helper function for findproc -- recursively determines k-size subsets of B[1:i] that sum to T
function f(i, T, k, memo = Dict())
    # Base cases
    if i >= 1 && T == 0 && k == 0
        Set()
    elseif T == 0 && k >= 1
        Nothing
    elseif i == 0 && T >= 1
        Nothing
    elseif T >= 1 && k == 0
        Nothing
    elseif T <= -1
        Nothing
    elseif i == 1
        T == k * B[1] ? Set([repeat([B[1]], k)]) : Nothing
    elseif i >= 1 && k == 1
        index = findall(x -> x == T, B[1:i])
        length(index) > 0 ? Set([[B[index[1]]]]) : Nothing
        # Recursive case with memoization
    else
        if haskey(memo, (i, T, k))
            get(memo, (i, T, k), 0)
        else
            set = unionF(f(i - 1, T, k, memo), mapCat(f(i, T - B[i], k - 1, memo), [B[i]]))
            get!(memo, (i, T, k), set)
        end
    end
end

# Helper function for f -- set unions with Nothing
function unionF(x, y)
    if x == Nothing && y == Nothing
        Nothing
    elseif x == Nothing
        y
    elseif y == Nothing
        x
    else
        union(x, y)
    end
end

# Helper function for f -- concatenate e with each element in set S
function mapCat(S, e)
    S == Nothing ? Nothing : Set([vcat(s, e) for s in S])
end

end


text:
"""
    ModelFactory

A collection of functions to generate the feedback networks and baseline models
used in evaluations and examples.
"""
module ModelFactory
using Flux
using ..Splitters
using ..Mergers
using ..FeedbackChains

export basenetparams, featnetparams, kernnetparams, slimnetparams, largenetparams
# export forwardconvnet, feedbackconvnet
export spoerer_model_fw, spoerer_model_b, spoerer_model_bf, spoerer_model_bk,
       spoerer_model_bl, spoerer_model_bt, spoerer_model_blt
################################################################################
# Parameter definitions
################################################################################
# baseline architectures from Spoerer et al. (2017)
basenetparams = (((3, 3), 1=>32, relu), ((3, 3), 32=>32, relu))
featnetparams = (((3, 3), 1=>45, relu), ((3, 3), 45=>45, relu))
kernnetparams = (((5, 5), 1=>32, relu), ((5, 5), 32=>32, relu))
# additional baseline architectures
slimnetparams = (((3, 3), 1=>17, relu), ((3, 3), 17=>17, relu), ((3, 3), 17=>17, relu),
                 ((3, 3), 17=>17, relu), ((3, 3), 17=>17, relu), ((3, 3), 17=>17, relu),
                 ((3, 3), 17=>17, relu), ((3, 3), 17=>17, relu))
largenetparams = (((3, 3), 1=>32, relu), ((3, 3), 32=>32, relu), ((3, 3), 32=>32, relu),
                  ((3, 3), 32=>32, relu), ((3, 3), 32=>32, relu), ((3, 3), 32=>32, relu),
                  ((3, 3), 32=>32, relu), ((3, 3), 32=>32, relu))

################################################################################
# Network generators
################################################################################
"""
    forwardconvnet(layerparams; classes=10, poolsize=(28, 28))

Construct a feedforward ConvNet like the ones used as baselines in Jarvers & Neumann (2019).
"""
function forwardconvnet(layerparams; classes=10, poolsize=(28, 28))
    # convolution layers according to parametrization
    layers = map(params -> Conv(params..., pad=params[1][1] .÷ 2), layerparams)
    # readout layer with global max-pooling and number of
    features = size(layers[end].weight)[4]
    poolglobal = MaxPool(poolsize)
    flatten(x) = reshape(x, :, size(x, 4))
    readout = Dense(features, classes)
    # build network
    return Chain(layers..., poolglobal, flatten, readout, softmax)
end # function forwardconvnet

"""
    feedbackconvnet(forwardparams, backwardparams, combineop; features=32,
        classes=10, poolsize=(28, 28), examplestate=zeros(Float32, 1,1,1,1))

Construct a ConvNet with feedback connections.
"""
function feedbackconvnet(forwardparams, backwardparams, combineop; features=32,
        classes=10, poolsize=(28, 28), examplestate=zeros(Float32, 1,1,1,1))
    layers = map(params -> Conv(params..., pad=params[1][1] .÷ 2), forwardparams)
    feedbacks = map(params -> Conv(params..., pad=params[1][1] .÷ 2), backwardparams)
    splits = map(i -> Splitter("split_$i"), 1:length(forwardparams))
    #merges = map((i, fb) -> )
    # TODO: complete this function
end # function feedbackconvnet

"""
    spoerer_model_fw(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)

Generate one of the forward models (B, B-K, B-F) from the paper:
  Spoerer, C.J., McClure, P. & Kriegeskorte, N. (2017).
  Recurrent convolutional neural networks: a better model of biological object recognition.
  Frontiers in Psychology 8, 1551.
"""
function spoerer_model_fw(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)
    return Chain(
        Conv(kernel, channels=>features, relu, pad=map(x -> x ÷ 2, kernel)),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool((2,2), stride=(2,2)),
        Conv(kernel, features=>features, relu, pad=map(x -> x ÷ 2, kernel)),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool(map(x -> x ÷ 2, inputsize)),
        flatten,
        Dense(features, classes, σ)
    )
end # function spoerer_model_fw

spoerer_model_b(T; channels=1, inputsize=(28,28), classes=10) =
    spoerer_model_fw(T, channels=channels, inputsize=inputsize, classes=classes)

spoerer_model_bk(T; channels=1, inputsize=(28,28), classes=10) =
    spoerer_model_fw(T, channels=channels, inputsize=inputsize, classes=classes, kernel=(5,5))

spoerer_model_bf(T; channels=1, inputsize=(28,28), classes=10) =
    spoerer_model_fw(T, channels=channels, inputsize=inputsize, classes=classes, features=64)

# TODO: currently, the forward, lateral and backward convolutions all have their
#       own biases. This should not make the model more powerful, as they are
#       combined additively before the non-linearity, but it wastes some resources.
"""
    spoerer_model_bl(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)

Generate the convolutional neural network with lateral recurrence (BL) from:
  Spoerer, C.J., McClure, P. & Kriegeskorte, N. (2017).
  Recurrent convolutional neural networks: a better model of biological object recognition.
  Frontiers in Psychology 8, 1551.
"""
function spoerer_model_bl(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)
    return FeedbackChain(
        Conv(kernel, channels=>features, pad=map(x -> x ÷ 2, kernel)),
        Merger("l1", ConvTranspose(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)), +),
        x -> relu.(x),
        Splitter("l1"),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool((2,2), stride=(2,2)),
        Conv(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)),
        Merger("l2", ConvTranspose(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)), +),
        x -> relu.(x),
        Splitter("l2"),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool(map(x -> x ÷ 2, inputsize)),
        flatten,
        Dense(features, classes, σ)
    )
end # function spoerer_model_bl

"""
    spoerer_model_bt(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)

Generate the convolutional neural network with top-down recurrence (BT) from:
  Spoerer, C.J., McClure, P. & Kriegeskorte, N. (2017).
  Recurrent convolutional neural networks: a better model of biological object recognition.
  Frontiers in Psychology 8, 1551.
"""
function spoerer_model_bt(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)
    return FeedbackChain(
        Conv(kernel, channels=>features, pad=map(x -> x ÷ 2, kernel)),
        Merger("l2",  ConvTranspose((2,2), features=>features, stride=2), +),
        x -> relu.(x),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool((2,2), stride=(2,2)),
        Conv(kernel, features=>features, relu, pad=map(x -> x ÷ 2, kernel)),
        Splitter("l2"),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool(map(x -> x ÷ 2, inputsize)),
        flatten,
        Dense(features, classes, σ)
    )
end # function spoerer_model_bt

"""
    spoerer_model_blt(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)

Generate the convolutional neural network with lateral and top-down recurrence (BLT) from:
  Spoerer, C.J., McClure, P. & Kriegeskorte, N. (2017).
  Recurrent convolutional neural networks: a better model of biological object recognition.
  Frontiers in Psychology 8, 1551.
"""
function spoerer_model_blt(T; channels=1, inputsize=(28, 28), kernel=(3,3), features=32, classes=10)
    return FeedbackChain(
        Conv(kernel, channels=>features, pad=map(x -> x ÷ 2, kernel)),
        Merger("l1", ConvTranspose(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)), +),
        Merger("l2", ConvTranspose((2,2), features=>features, stride=2), +),
        x -> relu.(x),
        Splitter("l1"),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool((2,2), stride=(2,2)),
        Conv(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)),
        Merger("l2", ConvTranspose(kernel, features=>features, pad=map(x -> x ÷ 2, kernel)), +),
        x -> relu.(x),
        Splitter("l2"),
        LRN(T(1.0), T(1.0), T(0.5), 1),
        MaxPool(map(x -> x ÷ 2, inputsize)),
        flatten,
        Dense(features, classes, σ)
    )
end # function spoerer_model_blt

################################################################################
# Helper functions
################################################################################
# helper function to turn Arrays and similar into tuples
tuplefy = (xs...) -> xs

# helper function to turn high-dimensional arrays into 2D (feature x batch)
flatten(x) = reshape(x, :, size(x, 4))

"""
    LRN{T,I}

Local response normalization layer. Input `i` is processed according to
out(x,y,f,b) = x * [b + α * sum( i(x, y, f-k÷2:f+k÷2, b)^2 )]^(-β)

Todo: β is currently ignored (always set to 0.5)
"""
struct LRN{T,I}
    b::T
    α::T
    β::T
    k::I
end

"""
    (l::LRN)(i)

Applies a local response normalization layer according to:
out(x,y,f,b) = x * [c + α * sum( i(x, y, f-k÷2:f+k÷2, b)^2 )]^(-β)

Todo: β is currently ignored (always set to 0.5)
"""
function (l::LRN)(x)
    ω = similar(x)
    fsize = size(x, 3)
    depth = l.k ÷ 2
    buffer = similar(x, size(x, 1), size(x, 2), 2*depth+1, size(x, 4))
    for i ∈ 1:fsize
        f_min = max(1, i - depth)
        f_max = min(fsize, i + depth)
        buffer = Flux.Tracker.data(x[:, :, f_min:f_max, :])
        ω[:, :, i, :] = sum(buffer.^2, dims=3)
    end
    return x ./ sqrt.(l.b .+ l.α .* ω)
end # function (l::LRN)

Flux.@treelike LRN

end # module ModelFactory


text:
using PyPlot
# plt[:style][:use]("classic")
# plt[:rcParams]["font.family"][1] = "serif"
# include("../src/DSR.jl")
using DSR

function test1()
    # A classic Single-Degree-of-Freedom system under resonant harmonic load
    T, zta = 1.0, 0.05
    omg = 2.0*pi/T
    t1, dt = 20.0, 0.005
    s = SDF(omg, zta)
    p = (t) -> sin(2.0*π/T*t)
    r = response(s,p,t1,dt)

    figure("Resonant Response",(9,3))
    plot(r[:,1], r[:,2])
    xlabel("Time [s]")
    ylabel("Displacement [m]")
    grid(true)
    # savefig("../Docs/resonant_response.png"; dpi=96, bbox_inches="tight")
    show()
end

function test2()
    # A classic Single-Degree-of-Freedom system under resonant harmonic load
    T, zta = 1.0, 0.05
    omg = 2.0*pi/T
    xi, alpha = 0.015, 0.2
    s = SDF_VD(omg, zta, xi, alpha)

    a_g = readtxt("EQ-S-1.txt")
    dt = 0.005
    nstep = 9000
    r  = response(s,  a_g, dt, dt, nstep)

    figure("Seismic Response of SDF with VD",(9,9))
    subplot(3,1,1)
    plot(r[:,1], r[:,2])
    xlabel("\$ t [\\rm{s}] \$")
    ylabel("\$ u [\\rm{m}] \$")
    grid(true)

    subplot2grid((3,1),(1,0),rowspan=2,colspan=1)
    plot(r[:,2],r[:,5],label="Inherent Damping")
    plot(r[:,2],r[:,6],label="Nonlinear Damper")
    xlabel("\$ u [\\rm{m}] \$")
    ylabel("\$ f_d/m [ \\rm{m/s^{2}} ] \$")
    grid(true)
    legend()
    # savefig("../Docs/seismic_response_VD.png"; dpi=96, bbox_inches="tight")
    show()
end

test1();
# test2();

max_stars_repo_name: Panchatantra/DSR

max_issues_repo_name: Panchatantra/DSR

max_forks_repo_name: Panchatantra/DSR


text:
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
# ==============================================================================
#
#   Tests related to the string cropping.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #

@testset "Left cropping" begin
    str = "Test 😅 \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"

    expected = "😅 \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"
    cropped_str = *(left_crop(str, 5)...)
    @test cropped_str == expected

    expected = "  \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"
    cropped_str = *(left_crop(str, 6)...)
    @test cropped_str == expected

    expected = "\e[38;5;231;48;5;243mst 😅 \e[38;5;201;48;5;243mTest\e[0m"
    cropped_str = *(left_crop(str, 10)...)
    @test cropped_str == expected

    expected = "\e[38;5;231;48;5;243m\e[38;5;201;48;5;243mest\e[0m"
    cropped_str = *(left_crop(str, 17)...)
    @test cropped_str == expected

    ansi_escape_seq, ~ = left_crop(str, 17)
    @test ansi_escape_seq == "\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m"
end

@testset "Fit field" begin
    str = "Test 😅 \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"
    printable_string_width = printable_textwidth(str)

    # Cropping from the right
    # ==========================================================================

    cropped_str = fit_string_in_field(str, 8)
    expected = "Test 😅…\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    @test cropped_str == expected

    cropped_str = fit_string_in_field(
        str,
        8;
        printable_string_width
    )
    expected = "Test 😅…\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    @test cropped_str == expected

    cropped_str = fit_string_in_field(str, 8; add_space_in_continuation_char = true)
    expected = "Test   …\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    @test cropped_str == expected

    cropped_str = fit_string_in_field(
        str,
        8;
        add_space_in_continuation_char = true,
        keep_ansi = false
    )
    expected = "Test   …"
    @test cropped_str == expected

    # Cropping from the left
    # ==========================================================================

    cropped_str = fit_string_in_field(str, 8; crop_side = :left)
    expected = "\e[38;5;231;48;5;243m…😅 \e[38;5;201;48;5;243mTest\e[0m"
    @test cropped_str == expected

    cropped_str = fit_string_in_field(
        str,
        8;
        crop_side = :left,
        printable_string_width
    )
    expected = "\e[38;5;231;48;5;243m…😅 \e[38;5;201;48;5;243mTest\e[0m"
    @test cropped_str == expected

    cropped_str = fit_string_in_field(
        str,
        8;
        add_space_in_continuation_char = true,
        crop_side = :left
    )
    expected = "\e[38;5;231;48;5;243m…   \e[38;5;201;48;5;243mTest\e[0m"
    @test cropped_str == expected
end

@testset "Right cropping" begin
    str = "Test 😅 \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"

    expected = "Test 😅 \e[38;5;231;48;5;243mTest 😅\e[38;5;201;48;5;243m\e[0m"
    cropped_str = *(right_crop(str, 5)...)
    @test cropped_str == expected

    expected = "Test 😅 \e[38;5;231;48;5;243mTest  \e[38;5;201;48;5;243m\e[0m"
    cropped_str = *(right_crop(str, 6)...)
    @test cropped_str == expected

    expected = "Test 😅\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    cropped_str = *(right_crop(str, 13)...)
    @test cropped_str == expected

    expected = "Test  \e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    cropped_str = *(right_crop(str, 14)...)
    @test cropped_str == expected

    ~, ansi_escape_seq = right_crop(str, 14)
    @test ansi_escape_seq == "\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"

    # Test passing the string printable width.
    printable_string_width = printable_textwidth(str)
    expected = "Test 😅 \e[38;5;231;48;5;243mTest 😅\e[38;5;201;48;5;243m\e[0m"
    cropped_str = *(right_crop(str, 5; printable_string_width)...)
    @test cropped_str == expected
end

@testset "Corner cases" begin
    str = "Test 😅 \e[38;5;231;48;5;243mTest 😅 \e[38;5;201;48;5;243mTest\e[0m"
    cropped_str = fit_string_in_field(str, 25)
    @test cropped_str == str

    expected = "…\e[38;5;231;48;5;243m\e[38;5;201;48;5;243m\e[0m"
    cropped_str = fit_string_in_field(str, 0)
    @test cropped_str == expected

    expected = "…"
    cropped_str = fit_string_in_field(str, 0; keep_ansi = false)
    @test cropped_str == expected
end