text:
## logging with dates
# JL Documentation page rendering may not compatible with LogginExtras, so use timestamp instead.
timestamp() = Dates.format(now(), "yyyy-mm-dd HH:MM:SS ")

## UUID generation
const UUID4 = uuid4(UUIDs.MersenneTwister(39871))

## simple functions
do_nothing() = nothing
do_nothing(x) = nothing
do_nothing(x, y) = nothing

isok(x::Nothing) = true
isok(x::Bool) = x
isok(x::AbstractString) = occursin(r"y(es)?|ok?|t(rue)?|^1$"i, x)
isok(x) = true  # default is true

## parse default inputs/outputs
function throw_invalid_xxputs(any)
	throw(ErrorException("DataTypeError: Elements of inputs and outputs can only be one of `name::Union{String,Symbol}`, `name => default_value`, `name => data_type`, `name => default_value => data_type`, `name => data_type => default_value`. Invalid value: $any"))
end

"""
    parse_default(v)

Parsing `inputs` and `outputs` when creating `Program` objects.

Return `xxputs::Vector{String}, xxput_types::Vector{DataType}, default_xxputs::Vector`.

## Valid `v` element types

- `name::String`: no default value.

- `name::String => value`: set default value, except `value` is `nothing` (default value not set).

- `name::String => value_type::DataType`: no default value, but value type.

- `name::String => value => value_type::DataType`: set default value and value type.

- `name::String => value_type::DataType => value`: set default value and value type.
"""
function parse_default(v::Vector{String})
    xxput_types = DataType[Any for i = 1:length(v)]
    default_xxputs = Nothing[nothing for i = 1:length(v)]
    return v, xxput_types, default_xxputs
end

function parse_default(v::Vector)
    n = length(v)
    xxputs = Vector{String}(undef, n)
    xxput_types = Vector{DataType}(undef, n)
    default_xxputs = Vector{Any}(undef, n)
    for (i, ele) in enumerate(v)
        name, type, default = parse_default_element(ele)
        xxputs[i] = name
        xxput_types[i] = type
        default_xxputs[i] = default
    end
    return xxputs, xxput_types, default_xxputs
end

parse_default(s::String) = parse_default([s])
parse_default(s::AbstractString) = parse_default([str(s)])
parse_default(p::Pair) = parse_default([p])
parse_default(d::Dict) = parse_default([p for p in d])
parse_default(s::Symbol) = parse_default([str(s)])
parse_default(any) = throw_invalid_xxputs(any)

function parse_default_element(ele::String)
    return ele, Any, nothing
end
function parse_default_element(ele::AbstractString)
    return str(ele), Any, nothing
end
function parse_default_element(ele::Symbol)
    return str(ele), Any, nothing
end
function parse_default_element(ele::Pair)
	if ele.first isa AbstractString || ele.first isa Symbol
	    name = str(ele.first)
		data_type, value = parse_arg_info(ele.second)
		name, data_type, value
	else
		throw_invalid_xxputs(ele)
	end
end
parse_default_element(any) = throw_invalid_xxputs(any)

function parse_arg_info(data_type::DataType)
	data_type, nothing
end
function parse_arg_info(p::Pair)
	parse_arg_info_pair(p.first, p.second)
end
function parse_arg_info(value)
	Any, value
end
function parse_arg_info_pair(a::DataType, b::DataType)
	if a isa b
		b, a
	elseif b isa a
		a, b
	else
		throw(ErrorException("DataTypeError: $(a) and $(b) are not inclusive."))
	end
end
parse_arg_info_pair(a, b::DataType) = b, convert_data_type(a, b)
parse_arg_info_pair(a::DataType, b) = a, convert_data_type(b, a)

function convert_data_type(value, data_type::DataType)
	if isa(value, data_type)
		return value
	else
		try
			convert(data_type, value)
		catch
			throw(ErrorException("DataTypeError: cannot convert $value to $data_type."))
		end
	end
end

## convert inputs/outputs to Dict{String} in `run(p, inputs, outputs)`
function to_xxput_dict(p::Pair{String, V}) where V
	Dict(p.first => p.second)
end
function to_xxput_dict(p::Pair)
	Dict(str(p.first) => p.second)
end
function to_xxput_dict(v::Vector{V}) where V <: Pair
	res = Dict{String,Any}()
	for p in v
		res[str(p.first)] = p.second
	end
	res
end
function to_xxput_dict(d::Dict)
	res = Dict{String,Any}()
	for p in v
		res[str(p.first)] = p.second
	end
	res
end
to_xxput_dict(d::Dict{String}) = d
to_xxput_dict(any) = throw(ErrorException("DataTypeError: cannot run Program: cannot convert to Dict{String}: inputs, outputs, or returned value of infer_outputs(inputs)"))

## String/Cmd conversion
"""
    to_str(x) -> String
    str(x) -> String

Convert `x` to `String`.

- `x::Cmd`: remove backticks (return `string(x)[2:end-1]`).

- `x::Nothing`: return `""`.

- `x::Vector`: join elements with `"_"` as delim.

- `x::Any`: return `string(x)`.
"""
str(x::Cmd) = string(x)[2:end-1]
str(::Nothing) = ""
str(x::Vector) = join(str.(x), "_")
str(x::Regex) = x.pattern
str(any) = string(any)
to_str = str

"""
    to_cmd(x) -> Cmd

Convert `x` to Cmd.

Exception: when `x::Nothing`, return `nothing::Nothing`.
"""
to_cmd(c::Cmd) = c
to_cmd(s::AbstractString) = Cmd([str(s)])
to_cmd(v::Vector{String}) = Cmd(v)
to_cmd(v::Vector{T}) where T = Cmd(str.(v))
to_cmd(::Nothing) = nothing
to_cmd(i::Number) = Cmd([str(i)])
to_cmd(any) = Cmd([str(any)])

ValidInputTypes = Union{Cmd, AbstractString, Vector, Number}

"""
    split(c::Cmd) = c.exec :: Vector{String}

Return splitted arguments of `Cmd`.
"""
function Base.split(c::Cmd)
    c.exec
end

## path functions
"""
    replaceext(path, replacement::AbstractString)

If the last component of a path contains a dot, leave everything before the dot as usual, and everything after the dot is replaced by `replacement`. Otherwise, `replacement` will be appended to `path`.

If `replacement` is empty, the last dot will be removed.
"""
function replaceext(path::String, replacement::AbstractString)
    a, b = splitext(path)
    if isempty(replacement)
        a
    elseif replacement[1] == '.'
        a * replacement
    else
        a * "." * replacement
    end
end
replaceext(path, replacement::AbstractString) = replaceext(to_str(path), replacement::String)

"""
    removeext(path)

If the last component of a path contains a dot, leave everything before the dot as usual, and everything including and after the dot is discarded.
"""
removeext(path::String) = splitext(path)[1]
removeext(path) = removeext(to_str(path))

## redirecting

# save the default std
stdout_origin = Base.stdout
stderr_origin = Base.stderr

"""
	set_default_stdout()
	set_default_stderr()

Backup the current stdout and stderr as the default stdout and stderr when redirecting stdout/stderr fails when calling `redirect_to_files`. It happens when an old stream is closed and then redirected to.

Please do not close them if you will use any of `run(::Program, ...)`, `redirect_to_files(...)`, `JobScheduler.Job(::Program, ...)`.

!!! note
    Redirecting in Julia are not thread safe, so unexpected redirection might be happen if you are running programs in different `Tasks` or multi-thread mode.
"""
function set_default_stdout()
	global stdout_origin = Base.stdout
end
function set_default_stderr()
	global stderr_origin = Base.stderr
end

# extend IO
Base.close(::Nothing) = nothing
Base.SimpleLogger(::Nothing) = nothing
Base.with_logger(f::Function, ::Nothing) = f()

Base.redirect_stdout(f::Function, ::Nothing) = f()
Base.redirect_stderr(f::Function, ::Nothing) = f()
Base.redirect_stdout(f::Function, ::Nothing) = f()

handle_open(::Nothing, mode) = nothing
handle_open(io::IO, mode) = io # do not change and do not close when exit
handle_open(file::AbstractString, mode) = open(file::AbstractString, mode)

"""
	redirect_to_files(f::Function, outfile; mode="a+")
	redirect_to_files(f::Function, outfile, errfile; mode="a+")
	redirect_to_files(f::Function, outfile, errfile, logfile; mode="a+")

Redirect outputs of function `f` to file(s).

- `xxxfile`: File path (`AbstractString`), `nothing` or `::IO`. `nothing` means no redirect. Files can be the same.
- `mode`: same as `open(..., mode)`.

Caution: If `xxxfile` is an `IO`, it won't be closed. Please use `close(io)` or `JobSchedulers.close_in_future(io, jobs)` manually!

!!! note
    Redirecting in Julia are not thread safe, so unexpected redirection might be happen if you are running programs in different `Tasks` or multi-thread mode.
"""
function redirect_to_files(f::Function, outfile, errfile, logfile; mode="a+")
	out = handle_open(outfile, mode)
	err = errfile == outfile ? out : handle_open(errfile, mode)
	log = logfile == outfile ? out : logfile == errfile ? err : handle_open(logfile, mode)

	old_stdout = Base.stdout
	old_stderr = Base.stderr

	isnothing(out) || redirect_stdout(out)
	isnothing(err) || redirect_stderr(err)

	# it seems with_logger is thread safe?
	logger = SimpleLogger(log)
	res = with_logger(logger) do
		try
			f()
		catch e
			@error "" exception=(e, catch_backtrace())
			e
		end
	end

	# before switch to old stdxxx, check whether it is opened. It might be happen due to redirect_xxx is not thread safe!
	if !isnothing(out)
		if isopen(old_stdout)
			try
				# it is not atomic, so use try
				redirect_stdout(old_stdout)
			catch
				redirect_stdout(stdout_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stdout(stdout_origin)
		end
	end

	if !isnothing(err)
		if isopen(old_stderr)
			try
				# it is not atomic, so use try
				redirect_stderr(old_stderr)
			catch
				redirect_stderr(stderr_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stderr(stderr_origin)
		end
	end

	outfile isa IO || close(out)
	errfile isa IO || close(err)
	logfile isa IO || close(log)
	if res isa Exception
		throw(res)
	else
		res
	end
end

function redirect_to_files(f::Function, outfile, errfile; mode="a+")
	out = handle_open(outfile, mode)
	err = errfile == outfile ? out : handle_open(errfile, mode)

	old_stdout = Base.stdout
	old_stderr = Base.stderr

	isnothing(out) || redirect_stdout(out)
	isnothing(err) || redirect_stderr(err)

	# it seems with_logger is thread safe?
	logger = SimpleLogger(err)
	res = with_logger(logger) do
		try
			f()
		catch e
			@error "" exception=(e, catch_backtrace())
			e
		end
	end

	# before switch to old stdxxx, check whether it is opened. It might be happen due to redirect_xxx is not thread safe!
	if !isnothing(out)
		if isopen(old_stdout)
			try
				# it is not atomic, so use try
				redirect_stdout(old_stdout)
			catch
				redirect_stdout(stdout_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stdout(stdout_origin)
		end
	end

	if !isnothing(err)
		if isopen(old_stderr)
			try
				# it is not atomic, so use try
				redirect_stderr(old_stderr)
			catch
				redirect_stderr(stderr_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stderr(stderr_origin)
		end
	end

	outfile isa IO || close(out)
	errfile isa IO || close(err)
	if res isa Exception
		throw(res)
	else
		res
	end
end

function redirect_to_files(f::Function, redirectfile; mode="a+")
	out = handle_open(redirectfile, mode)

	old_stdout = Base.stdout
	old_stderr = Base.stderr

	if !isnothing(out)
		redirect_stdout(out)
		redirect_stderr(out)
	end

	logger = SimpleLogger(out)

	res = with_logger(logger) do
		try
			f()
		catch e
			@error "" exception=(e, catch_backtrace())
		end
	end

	# before switch to old stdxxx, check whether it is opened. It might be happen due to redirect_xxx is not thread safe!
	if !isnothing(out)
		if isopen(old_stdout)
			try
				# it is not atomic, so use try
				redirect_stdout(old_stdout)
			catch
				redirect_stdout(stdout_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stdout(stdout_origin)
		end

		if isopen(old_stderr)
			try
				# it is not atomic, so use try
				redirect_stderr(old_stderr)
			catch
				redirect_stderr(stderr_origin) # defined when using Pipelines, can define manually by using set_default_stdout()
			end
		else
			redirect_stderr(stderr_origin)
		end
	end

	redirectfile isa IO || close(out)
	if res isa Exception
		throw(res)
	else
		res
	end
end


text:
println("I am being built...")
#run(`julia docs/make.jl`)

ext: jl


text:
export rule

@rule typeof(dot)(:in2, Marginalisation) (m_out::UnivariateNormalDistributionsFamily, m_in1::PointMass{ <: AbstractVector }, meta::AbstractCorrection) = begin
    x  = mean(m_in1)
    xi = x * weightedmean(m_out)
    w  = correction!(meta, x * precision(m_out) * x')
    return MvGaussianWeightedMeanPrecision(xi, w)
end


text:
module How
    export howmuch

    function howmuch(m::Int, n::Int)::Array{Tuple{String,String,String}}
        min_value, max_value = minmax(m,n)
        # Only values "37+63idx" are possible
        idx_to_value(idx) = 37 + 63idx
        value_to_idx(v) = (v - 37) / 63
        min_idx = ceil(Int,value_to_idx(min_value))
        max_idx = floor(Int,value_to_idx(max_value))
        result = []
        for idx = min_idx:max_idx
            money = idx_to_value(idx)
            boat_cost = money ÷ 7
            car_cost = money ÷ 9
            push!(result, ("M: $money","B: $boat_cost","C: $car_cost") )
        end
        result
    end

end

max_stars_repo_name: UlrichBerntien/Codewars-Katas

max_issues_repo_name: UlrichBerntien/Codewars-Katas

max_forks_repo_name: UlrichBerntien/Codewars-Katas


text:
function FFT_prepare(f)
	for _i in 1:f.len
		p=0
		i=_i-1
		for j in 1:f.bits
			p=(p*2)+ ( (i%2)==1 ? 1 : 0 )
			i=fld(i,2)
		end
		#println(_i-1," ",p)
		f.exchangeMap[_i]=p+1
	end
	for i in 1:fld(f.len,2)
		ang=-2*pi*(i-1)/f.len
		f.e[i]=cos(ang)+sin(ang)*im
	end
	return f
end

struct FFT
	bits::Integer
	len::Integer
	exchangeMap
	e
	FFT(bits)=FFT_prepare(new(	bits,
					2^bits,
					zeros(Int,2^bits),
					zeros(Complex,2^(bits-1))))
end

function (fft::FFT)(f)
	F=zeros(Complex,2^(fft.bits) )
	for i in 1:fft.len
		F[fft.exchangeMap[i]]=f[i]
	end

	for i in 0:(fft.bits-1)
		unitSize_2=2^i
		unitSize=2*unitSize_2
		unitNum=2^(fft.bits-i-1)

		#=
		for j in 1:unitSize:
			coef=fft.e[(j-1)*unitNum+1]
			for k in 0:(unitNum-1):
				pPre= unitSize*k+j
				pPost=pPre+unitSize_2
				pre=F[pPre]
				post=F[pPost]
				F[pPre]=pre+coef*post
				F[pPost]=pre-coef*post
		=#
		for k in 0:(unitNum-1)
			offset=unitSize*k+1
			for j in 0:(unitSize_2-1)
				coef=fft.e[j*unitNum+1]
				pPre=offset+j
				pPost=pPre+unitSize_2
				pre=F[pPre]
				post=F[pPost]*coef
				#println([pre,post,coef])
				F[pPre]=pre+post
				F[pPost]=pre-post
			end
		end
	end
	return F./fft.len
end

f2=FFT(2)
f2=FFT(2)
println(f2)
println(f2([0,1,2,3]))

println()
f3=FFT(3)
println(f3)
println(f3(0:7))
F10=FFT(10)
for i in 1:(10)
	F10(rand(1024))
end


text:
struct RODECompositeSolution{T,N,uType,uType2,EType,tType,randType,P,A,IType} <: DiffEqBase.AbstractRODESolution{T,N}
  u::uType
  u_analytic::uType2
  errors::EType
  t::tType
  W::randType
  prob::P
  alg::A
  interp::IType
  alg_choice::Vector{Int}
  dense::Bool
  tslocation::Int
  retcode::Symbol
  seed::UInt64
end
(sol::RODECompositeSolution)(t,deriv::Type=Val{0};idxs=nothing,continuity=:left) = sol.interp(t,idxs,deriv,sol.prob.p,continuity)
(sol::RODECompositeSolution)(v,t,deriv::Type=Val{0};idxs=nothing,continuity=:left) = sol.interp(v,t,idxs,deriv,sol.prob.p,continuity)

function DiffEqBase.build_solution(
        prob::DiffEqBase.AbstractRODEProblem,
        alg::Union{StochasticDiffEqCompositeAlgorithm,
                   StochasticDiffEqRODECompositeAlgorithm},t,u;
        W=[],timeseries_errors=length(u)>2,
        dense=false,dense_errors=dense,
        calculate_error = true,
        interp = LinearInterpolation(t,u),
        alg_choice=[1], retcode = :Default,seed = UInt64(0),
        kwargs...)

  T = eltype(eltype(u))
  if typeof(prob.u0) <: Tuple
    N = length((size(ArrayPartition(prob.u0))..., length(u)))
  else
    N = length((size(prob.u0)..., length(u)))
  end

  if typeof(prob.f) <: Tuple
    f = prob.f[1]
  else
    f = prob.f
  end

  if DiffEqBase.has_analytic(f)
    u_analytic = Vector{typeof(prob.u0)}()
    errors = Dict{Symbol,real(eltype(prob.u0))}()
    sol = RODECompositeSolution{T,N,typeof(u),typeof(u_analytic),typeof(errors),typeof(t),typeof(W),
                       typeof(prob),typeof(alg),typeof(interp)}(u,u_analytic,errors,t,W,prob,alg,interp,alg_choice,dense,0,retcode,seed)
    if calculate_error
      calculate_solution_errors!(sol;timeseries_errors=timeseries_errors,dense_errors=dense_errors)
    end
    return sol
  else
    return RODECompositeSolution{T,N,typeof(u),typeof(nothing),typeof(nothing),typeof(t),typeof(W),
                       typeof(prob),typeof(alg),typeof(interp)}(u,nothing,nothing,t,W,prob,alg,interp,alg_choice,dense,0,retcode,seed)
  end
end

function DiffEqBase.solution_new_retcode(sol::RODECompositeSolution{T,N,uType,uType2,EType,tType,randType,P,A,IType},retcode) where {T,N,uType,uType2,EType,tType,randType,P,A,IType}
  RODECompositeSolution{T,N,uType,uType2,EType,tType,randType,P,A,IType}(
                        sol.u,sol.u_analytic,sol.errors,sol.t,sol.W,sol.prob,
                        sol.alg,sol.interp,sol.alg_choice,sol.dense,sol.tslocation,
                        retcode,sol.seed)
end

function DiffEqBase.solution_new_tslocation(sol::RODECompositeSolution{T,N,uType,uType2,EType,tType,randType,P,A,IType},tslocation) where {T,N,uType,uType2,EType,tType,randType,P,A,IType}
  RODECompositeSolution{T,N,uType,uType2,EType,tType,randType,P,A,IType}(
                        sol.u,sol.u_analytic,sol.errors,sol.t,sol.W,sol.prob,
                        sol.alg,sol.interp,sol.alg_choice,sol.dense,tslocation,
                        sol.retcode,sol.seed)
end


text:
# Use
#    @warnpcfail precompile(args...)
# if you want to be warned when a precompile directive fails
macro warnpcfail(ex::Expr)
    modl = __module__
    file = __source__.file === nothing ? "?" : String(__source__.file)
    line = __source__.line
    quote
        $(esc(ex)) || @warn """precompile directive
     $($(Expr(:quote, ex)))
 failed. Please report an issue in $($modl) (after checking for duplicates) or remove this directive.""" _file=$file _line=$line
    end
end

const __bodyfunction__ = Dict{Method,Any}()

# Find keyword "body functions" (the function that contains the body
# as written by the developer, called after all missing keyword-arguments
# have been assigned values), in a manner that doesn't depend on
# gensymmed names.
# `mnokw` is the method that gets called when you invoke it without
# supplying any keywords.
function __lookup_kwbody__(mnokw::Method)
    function getsym(arg)
        isa(arg, Symbol) && return arg
        @assert isa(arg, GlobalRef)
        return arg.name
    end

    f = get(__bodyfunction__, mnokw, nothing)
    if f === nothing
        fmod = mnokw.module
        # The lowered code for `mnokw` should look like
        #   %1 = mkw(kwvalues..., #self#, args...)
        #        return %1
        # where `mkw` is the name of the "active" keyword body-function.
        ast = Base.uncompressed_ast(mnokw)
        if isa(ast, Core.CodeInfo) && length(ast.code) >= 2
            callexpr = ast.code[end-1]
            if isa(callexpr, Expr) && callexpr.head == :call
                fsym = callexpr.args[1]
                if isa(fsym, Symbol)
                    f = getfield(fmod, fsym)
                elseif isa(fsym, GlobalRef)
                    if fsym.mod === Core && fsym.name === :_apply
                        f = getfield(mnokw.module, getsym(callexpr.args[2]))
                    elseif fsym.mod === Core && fsym.name === :_apply_iterate
                        f = getfield(mnokw.module, getsym(callexpr.args[3]))
                    else
                        f = getfield(fsym.mod, fsym.name)
                    end
                else
                    f = missing
                end
            else
                f = missing
            end
        else
            f = missing
        end
        __bodyfunction__[mnokw] = f
    end
    return f
end

function _precompile_()
    ccall(:jl_generating_output, Cint, ()) == 1 || return nothing
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:name, :study_id),Tuple{Nothing,Nothing}},Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},FanBeamGeometry{Float32,DefaultTomograph},CTImage{Float32,Array{Float32,2}}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:name, :study_id),Tuple{Nothing,Nothing}},Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},FanBeamGeometry{Float64,DefaultTomograph},CTImage{Float64,Array{Float64,2}}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:name, :study_id),Tuple{Nothing,Nothing}},Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},ParallelBeamGeometry{Float32,DefaultTomograph},CTImage{Float32,Array{Float32,2}}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:name, :study_id),Tuple{Nothing,Nothing}},Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},ParallelBeamGeometry{Float64,DefaultTomograph},CTImage{Float64,Array{Float64,2}}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:name, :study_id),Tuple{Symbol,String}},Type{CTScanner},FBPInfo{ParallelBeamGeometry{Float32,DefaultTomograph},FBP{MartaCT.RadonAlgorithm.Filters.RamLak}},MartaCT.CTScan.CTImageData.CTData{Float32,Array{Float32,2}}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:nϕ,),Tuple{Int64}},Type{FanBeamGeometry},CircleImage{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:nϕ,),Tuple{Int64}},Type{FanBeamGeometry},CircleImage{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:nϕ,),Tuple{Int64}},Type{ParallelBeamGeometry},CircleImage{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:nϕ,),Tuple{Int64}},Type{ParallelBeamGeometry},CircleImage{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{CircleImage},Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{CircleImage},Type{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{GrayScaleLine},Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{GrayScaleLine},Type{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{GrayScalePyramid},Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{GrayScalePyramid},Type{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{WhiteRect},Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(Type)),NamedTuple{(:width,),Tuple{Int64}},Type{WhiteRect},Type{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(circle_image)),NamedTuple{(:radius, :calibration_value, :background, :rows),Tuple{Int64,Float32,Float32,Int64}},typeof(circle_image),Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(circle_image)),NamedTuple{(:radius, :calibration_value, :background, :rows),Tuple{Int64,Float64,Float64,Int64}},typeof(circle_image),Type{Float64}})
    Base.precompile(Tuple{Core.kwftype(typeof(gray_scale_image)),NamedTuple{(:swidth, :sheight, :gray_scale, :background),Tuple{Int64,Int64,IntervalSets.Interval{:closed,:closed,Float32},Float32}},typeof(gray_scale_image),Type{Float32}})
    Base.precompile(Tuple{Core.kwftype(typeof(gray_scale_image)),NamedTuple{(:swidth, :sheight, :gray_scale, :background),Tuple{Int64,Int64,IntervalSets.Interval{:closed,:closed,Float64},Float64}},typeof(gray_scale_image),Type{Float64}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},CircleImage{Float64}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},GrayScaleLine{Float64}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},GrayScalePyramid{Float32}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},GrayScalePyramid{Float64}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},WhiteRect{Float32}})
    Base.precompile(Tuple{Type{CTScanner{FBPInfo{G,A} where A<:AbstractFBP,D} where D where G},WhiteRect{Float64}})
    Base.precompile(Tuple{typeof(MartaCT.Utils._wrap_angle),Float32,Int64})
    Base.precompile(Tuple{typeof(MartaCT.Utils._wrap_angle),Float64,Int64})
    Base.precompile(Tuple{typeof(combine_images),ImageParams{Float32},Array{Float32,2},Array{Float32,2}})
    Base.precompile(Tuple{typeof(combine_images),ImageParams{Float64},Array{Float64,2},Array{Float64,2}})
    Base.precompile(Tuple{typeof(project_image),Some{CTImage{Float32,Array{Float32,2}}},RadonInfo{FanBeamGeometry{Float32,DefaultTomograph},Radon}})
    Base.precompile(Tuple{typeof(project_image),Some{CTImage{Float32,Array{Float32,2}}},RadonInfo{ParallelBeamGeometry{Float32,DefaultTomograph},Radon}})
    Base.precompile(Tuple{typeof(project_image),Some{CTImage{Float64,Array{Float64,2}}},RadonInfo{FanBeamGeometry{Float64,DefaultTomograph},Radon}})
    Base.precompile(Tuple{typeof(project_image),Some{CTImage{Float64,Array{Float64,2}}},RadonInfo{ParallelBeamGeometry{Float64,DefaultTomograph},Radon}})
    Base.precompile(Tuple{typeof(similar),MartaCT.CTScan.CTImageData.CTData{Float32,Array{Float32,2}},CTSinogram{Float32,Array{Float32,2}}})
    Base.precompile(Tuple{typeof(similar),MartaCT.CTScan.CTImageData.CTData{Float64,Array{Float64,2}},CTSinogram{Float64,Array{Float64,2}}})
    let fbody = try __lookup_kwbody__(which(polar2cart, (Array{Float32,2},StepRangeLen{Float32,Float64,Float64},StepRangeLen{Float32,Float64,Float64},))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Int64,Bool,Float32,Bool,MartaCT.Interpolation.NoInterpolation,typeof(polar2cart),Array{Float32,2},StepRangeLen{Float32,Float64,Float64},StepRangeLen{Float32,Float64,Float64},))
        end
    end
    let fbody = try __lookup_kwbody__(which(polar2cart, (Array{Float32,2},StepRangeLen{Float32,Float64,Float64},StepRangeLen{Float32,Float64,Float64},))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Int64,Bool,Float32,Bool,Nothing,typeof(polar2cart),Array{Float32,2},StepRangeLen{Float32,Float64,Float64},StepRangeLen{Float32,Float64,Float64},))
        end
    end
    let fbody = try __lookup_kwbody__(which(polar2cart, (Array{Float64,2},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Int64,Bool,Float64,Bool,MartaCT.Interpolation.NoInterpolation,typeof(polar2cart),Array{Float64,2},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},))
        end
    end
    let fbody = try __lookup_kwbody__(which(polar2cart, (Array{Float64,2},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Int64,Bool,Float64,Bool,Nothing,typeof(polar2cart),Array{Float64,2},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},StepRangeLen{Float64,Base.TwicePrecision{Float64},Base.TwicePrecision{Float64}},))
        end
    end
end


text:
using Test
using SpecialFunctions
using PyPlot
using FastGaussQuadrature

include("onTheFlyGH.jl")

@testset "Gauss-Hermite on the fly" begin
    compare()
end

include("GaussLaguerrePlots.jl")


text:
function unsafe_getbyoffset(::Type{T}, s::S, offset::Int) where {T, S}
    rt = Ref{T}()
    rs = Ref{S}(s)
    GC.@preserve rt rs begin
        pt = Ptr{UInt8}(Base.unsafe_convert(Ref{T}, rt))
        ps = Ptr{UInt8}(Base.unsafe_convert(Ref{S}, rs)) + offset
        Base._memcpy!(pt, ps, sizeof(T))
    end
    return rt[]
end

function check_getbyoffset(T, s, offset)
    S = typeof(s)
    isbitstype(T) || throw(ArgumentError("Can only cast into bitstype."))
    isbitstype(S) || throw(ArgumentError("Can only cast from bitstype."))
    @assert offset >= 0
    sizeof(T) + offset <= sizeof(S) || throw(ArgumentError("""
            sizeof(T) + offset <= sizeof(typeof(s)) must hold. Got:
            T = $(T), sizeof(T) = $(sizeof(T))
            S = $(S), sizeof(S) = $(sizeof(S))
            offset = $offset
            """)
    )
end

function getbyoffset(T, s, offset)
    @boundscheck check_getbyoffset(T, s, offset)
    unsafe_getbyoffset(T, s, offset)
end

function unsafe_setbyoffset(t::T, s::S, offset::Int) where {T, S}
    rt = Ref{T}(t)
    rs = Ref{S}(s)
    GC.@preserve rt rs begin
        pt = Ptr{UInt8}(Base.unsafe_convert(Ref{T}, rt)) + offset
        ps = Ptr{UInt8}(Base.unsafe_convert(Ref{S}, rs))
        Base._memcpy!(pt, ps, sizeof(T))
    end
    return rt[]::T
end

function check_setbyoffset(t, s, offset)
    T = typeof(t)
    S = typeof(s)
    isbitstype(T) || throw(ArgumentError("Can only set into bitstype."))
    isbitstype(S) || throw(ArgumentError("Can only set from bitstype."))
    @assert offset >= 0
    sizeof(T) >= sizeof(S) + offset || throw(ArgumentError("""
            sizeof(T) >= sizeof(typeof(s)) + offset must hold. Got:
            T = $(T), sizeof(T) = $(sizeof(T))
            S = $(S), sizeof(S) = $(sizeof(S))
            offset = $offset
            """)
    )
end

function setbyoffset(t, s, offset)
    @boundscheck check_setbyoffset(t, s, offset)
    unsafe_setbyoffset(t, s, offset)
end


text:
using Test
using SumOfSquares
using DynamicPolynomials

function quartic_constant_test(optimizer,
                   config::MOIT.TestConfig,
                   cone::SumOfSquares.PolyJuMP.PolynomialSet)
    atol = config.atol
    rtol = config.rtol

    @polyvar x
    p = x^4 + 4

    model = _model(optimizer)
    @variable(model, γ)

    cref = @constraint(model, p - γ in cone, basis=FixedPolynomialBasis([x^2]))

    @objective(model, Max, γ)

    optimize!(model)

    @test termination_status(model) == MOI.OPTIMAL

    @test primal_status(model) == MOI.FEASIBLE_POINT

    @test objective_value(model) ≈ 4 atol=atol rtol=rtol

    p = gram_matrix(cref)
    @test p isa SumOfSquares.GramMatrix
    @test getmat(p) ≈ ones(1, 1) atol=atol rtol=rtol
    @test p.basis isa FixedPolynomialBasis
    @test p.basis.polynomials == [x^2]

    S = SumOfSquares.SOSPolynomialSet{
        SumOfSquares.FullSpace, Monomial{true}, MonomialVector{true},
        SumOfSquares.Certificate.FixedBasis{typeof(cone),typeof(p.basis)}
    }
    @test list_of_constraint_types(model) == [(Vector{AffExpr}, S)]
    test_delete_bridge(
        model, cref, 1,
        ((MOI.VectorOfVariables, MOI.Nonnegatives, 0),
         (MOI.VectorOfVariables, MOI.PositiveSemidefiniteConeTriangle, 0)
        ))
end
sos_quartic_constant_test(optimizer, config)   = quartic_constant_test(optimizer, config, SOSCone())
sd_tests["sos_quartic_constant"] = sos_quartic_constant_test
sdsos_quartic_constant_test(optimizer, config) = quartic_constant_test(optimizer, config, SDSOSCone())
soc_tests["sdsos_quartic_constant"] = sdsos_quartic_constant_test
dsos_quartic_constant_test(optimizer, config)  = quartic_constant_test(optimizer, config, DSOSCone())
linear_tests["dsos_quartic_constant"] = dsos_quartic_constant_test


text:
# source:
# https://www.mathworks.com/help/pdf_doc/otherdocs/ode_suite.pdf

struct Ros23 <: AbstractRosenbrockSolver
end

function order(::Ros23)::Int
    return 2
end

function n_intermediates(::Ros23)
    return 3
end

function n_temporaries(::Ros23)
    return 3 + 4
end

function stepper!(ode_solver::Ros23, ode_fcn, basic_options::BasicOdeOptions, stepper_temp)

    # renaming for better readability
    t          = stepper_temp.t_old;
	x          = stepper_temp.x_old;
	dt         = stepper_temp.dt;
    n_prob_dim = n_dim(stepper_temp)

    # constants
    d = 1 / (2 + sqrt(2))
    e32 = 6 + sqrt(2)

    # inherent variables:
    J    = stepper_temp.jac # zeros(n_dim, n_dim)
    dfdt = stepper_temp.dfdt

    # temporaries
    # maximally 4 temporaries are needed

    k1 = view( stepper_temp.der_x, 1:n_prob_dim, 1 )
    k2 = view( stepper_temp.der_x, 1:n_prob_dim, 2 )
    k3 = view( stepper_temp.der_x, 1:n_prob_dim, 3 )

    f0 = view( stepper_temp.tmp, 1:n_prob_dim, 1 )  # f0 = similar(x)
    f1 = view( stepper_temp.tmp, 1:n_prob_dim, 2 )  # f1 = similar(x)
    f2 = view( stepper_temp.tmp, 1:n_prob_dim, 3 )  # f2 = similar(x)

    f1_dt   = view( stepper_temp.tmp, 1:n_prob_dim, 4 )
    h       = view( stepper_temp.tmp, 1:n_prob_dim, 5 )
    dfdq1   = view( stepper_temp.tmp, 1:n_prob_dim, 6 )
    dfdq2   = view( stepper_temp.tmp, 1:n_prob_dim, 7 )
    x_step  = stepper_temp.x_tmp

    # reusable temporary, but is needed one down there
    rhs   = view( stepper_temp.tmp, 1:n_prob_dim, 7 ) # similar(x)
    x_tmp = stepper_temp.x_tmp # similar(x)

    copy!(x_tmp, x)
    copy!(x_step, x)

    # derivatives
    ode_fcn(f0, t, x)

    # this needs h, x_step, dfdq1, dfdq2 -> 4 temporaries
    # inherent: J
    J    = numjac!(basic_options, ode_fcn, t, x, h, x_step, dfdq1, dfdq2, J)
    # temporary: f1_dt
    # inherent: dfdt
    dfdt = numder!(basic_options, ode_fcn, t, dt, x, f0, f1_dt, dfdt)

    # lu decomposition
    J .*= -d * dt
    for ii = 1:n_prob_dim
        J[ii, ii] += 1.0
    end
    luc = lu(J)

    # stages
    @. rhs = f0 + dt * d * dfdt
    k1 .= luc \ rhs

    @. x_tmp = x + dt * 0.5 * k1
    t1      = t + 0.5 * dt
    ode_fcn(f1, t1, x_tmp)
    @. rhs = f1 - k1
    k2  .= (luc \ rhs) + k1

    # new step
    stepper_temp.x_new = x + dt * k2
    stepper_temp.t_new = t + dt

    ode_fcn(f2, stepper_temp.t_new, stepper_temp.x_new)
    @. rhs = f2 - e32 * (k2 - f1) - 2 * (k1 - f0) + dt * d * dfdt
    k3 .= luc \ rhs

    # error estimation
    stepper_temp.x_est = stepper_temp.x_new + dt / 6 * (k1 - 2 * k2 + k3)

end

function numder!(basic_options, ode_fcn, t, dt, x, f0, f1_dt, dfdt)
    # h = sqrt(eps(abs(t)))
    h = 0.1 * dt
    ode_fcn(f1_dt, t + h, x)
    @. dfdt = (f1_dt - f0) / h
    return dfdt
end

function numjac!(basic_options, ode_fcn, t, x,
    h, x_step, dfdq1, dfdq2, J)

    n_dim = length(x)
    h .= basic_options.rel_tol .* abs.(x) .+ basic_options.abs_tol

    for ii = 1:n_dim

        x_step[ii] += h[ii]/2
        ode_fcn(dfdq2, t, x_step)

        x_step[ii] -= h[ii]
        ode_fcn(dfdq1, t, x_step)

        x_step[ii] += h[ii]/2

        J[:,ii] = (dfdq2 - dfdq1) ./ h[ii]
    end

    return J
end

#=
function numjac!(ode_fcn, t, x, f0, J, x1, f1)

    n = length(x)
    hs = 10 .* sqrt.(eps.(abs.(x)))

    for ii = 1:n
        x1[ii] += hs[ii]
        ode_fcn(f1, t, x1)
        J[:,ii] = (f1 - f0) ./ hs[ii]
        x1[ii] -= hs[ii]
    end

    return J
end
=#

max_stars_repo_name: HomoModelicus/julia

max_issues_repo_name: HomoModelicus/julia

max_forks_repo_name: HomoModelicus/julia


text:
using LinearAlgebraicRepresentation
using SparseArrays
Lar = LinearAlgebraicRepresentation
L = Lar
using ViewerGL
GL = ViewerGL

function twocubegrids(n,m,p)
    V,(VV,EV,FV,CV) = Lar.cuboidGrid([n,m,p],true)
    mybox = (V,CV,FV,EV)

    twocubs = Lar.Struct([mybox, L.t(.3,.4,.5), L.r(pi/5,0,0), L.r(0,0,pi/12), mybox])

    V,CV,FV,EV = Lar.struct2lar(twocubs)
    GL.VIEW([ GL.GLGrid(V,FV, GL.COLORS[1]) ]);

    cop_EV = Lar.coboundary_0(EV::Lar.Cells);
    cop_EW = convert(Lar.ChainOp, cop_EV);
    cop_FE = Lar.coboundary_1(V, FV::Lar.Cells, EV::Lar.Cells);
    W = convert(Lar.Points, V');

    V, copEV, copFE, copCF = Lar.Arrangement.spatial_arrangement( W, cop_EW, cop_FE)

    EV = Lar.cop2lar(copEV)
    FE = [findnz(copFE[k,:])[1] for k=1:size(copFE,1)]
    FV = [collect(Set(cat(EV[e] for e in FE[f]))) for f=1:length(FE)]
    FV = convert(Lar.Cells, FV)
    W = convert(Lar.Points, V')
    WW = [[k] for k=1:size(W,2)]

    #GL.VIEW(GL.numbering(.2)((W,[WW,EV])));

    V,CVs,FVs,EVs = Lar.pols2tria(W, copEV, copFE, copCF)
    GL.VIEW(GL.GLExplode(V,FVs,1.5,1.5,1.5,99,1));
    GL.VIEW(GL.GLExplode(V,EVs,1.5,1.5,1.5,99,1));
    GL.VIEW(GL.GLExplode(V,CVs,5,5,5,99,1));

    return V, copEV, copFE, copCF
end

twocubegrids(5,5,5);


text:
################################################################################

module CMD

using ..Die

abstract type Command end
abstract type AbortingCommand <: Command end
abstract type InputCommand <: Command end

struct ResetPhase <: AbortingCommand end
struct Exit <: AbortingCommand end
struct Phase <: AbortingCommand
	nbr
	max_nbr
end
Phase(max_nbr) = Phase(nothing,max_nbr)

struct Repeat <: Command end
struct Undo <: Command end
struct Help <: Command end

struct True <: InputCommand end
struct False <: InputCommand end
struct Option <: InputCommand opt end
struct Blank <: InputCommand end
struct Dice <: InputCommand dice end
Dice() = Dice(nothing)

########################################
Base.string(c::Command) = split(lowercase(string(typeof(c))), '.')[end]
Base.string(o::Option) = string(o.opt)

function Base.:(==)(a::Dice,b::Dice)
	isnothing(a.dice) && return true
	isnothing(b.dice) && return true
	return a.dice == b.dice
end

########################################
minmatch(::Command) = 1
minmatch(::AbortingCommand) = 0
minmatch(::ResetPhase) = 5
minmatch(::Option) = 0

function parse(cmd::Command, s::AbstractString)
	str = string(cmd)
	length(s) > length(str) && return nothing
	matchlen = minmatch(cmd) == 0 ? length(str) : max(length(s), minmatch(cmd))
	s == str[1:matchlen] ? cmd : nothing
end

parse(cmd::Dice, s::AbstractString) = (d = Die.parse(s); isnothing(d) ? nothing : Dice(d))
parse(cmd::Blank, s::AbstractString) = (isempty(strip(s)) ? Blank() : nothing)

function parse(cmd::Phase, s::AbstractString)
	s_parts = split(s)
	length(s_parts) != 2 && return nothing
	s_parts[1] != "phase" && return nothing

	nbr = tryparse(Int, s_parts[2])
	isnothing(nbr) && return nothing
	!(1 <= nbr <= cmd.max_nbr) && return nothing

	return Phase(nbr, cmd.max_nbr)
end

########################################
function parse(cmds::Vector{T}, s::AbstractString) where {T <: Command}
	s = lowercase(s)

	for cmd in cmds
		match = parse(cmd, s)
		!isnothing(match) && return match
	end
	return nothing
end

end

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

struct IOParser
	cmds::Vector{CMD.Command}
	in::IO
	out::IO

	IOParser(cmds) = new(cmds, stdin, stdout)
end

Base.print(iop::IOParser,a...) = print(iop.out, a...)
Base.println(iop::IOParser,a...) = println(iop.out, a...)
function Base.readline(iop::IOParser)
	str = readline(iop.in)
	!isopen(iop.in) && error("Input stream closed")
	return strip(str)
end

function read_input(iop::IOParser, options, prompt::String)
	print(iop, prompt)

	str = readline(iop)

	cmd = CMD.parse([iop.cmds; options], str)
	!isnothing(cmd) && return cmd

	!isempty(str) && println(iop,"Invalid input.")
	return read_input(iop, options, prompt)
end

function read_input(iop::IOParser, options::Vector{<:CMD.Command}, silent_options::Vector{<:CMD.Command})
	prompt = "["*strvec2str(string.(options), '/')*"] > "
	return read_input(iop, [options; silent_options], prompt)
end

function read_input(iop::IOParser, options::Vector{CMD.Dice}, silent_options::Vector{<:CMD.Command})
	prompt = "[$(strvec2str(Die.char.(instances(Die.Face)),','))] > "
	return read_input(iop, [options; silent_options], prompt)
end

function read_input(iop::IOParser, options::Vector{CMD.Blank}, silent_options::Vector{<:CMD.Command})
	return read_input(iop, [options; silent_options], "[Press enter to continue] > ")
end

function display_message(iop::IOParser,msg,header="-"^10)
	println(iop, header)
	msg = strip(msg)
	for p in split(msg, '\n')
		m = match(r"^(?:\s*[0-9]+\.|\s*-|And,|Or,|.+::)\s+", p)
		list_item_indent = isnothing(m) ? 0 : length(m.match)
		println(iop, wrap(p, width=60, subsequent_indent=' '^list_item_indent))
	end
end

max_stars_repo_name: mvmorin/queller-bot

max_issues_repo_name: mvmorin/queller-bot

max_forks_repo_name: mvmorin/queller-bot


text:
# This file is part of InvertedFiles.jl

using SimilaritySearch, LinearAlgebra, SparseArrays
using Base.Threads: SpinLock, @threads

export AbstractInvertedFile

"""
    abstract type AbstractInvertedFile <: AbstractSearchContext end

Abstract inverted file, actual data structures are [`WeightedInvertedFile`](@ref) and [`BinaryInvertedFile`](@ref)
"""
abstract type AbstractInvertedFile <: AbstractSearchContext end

"""
    length(idx::AbstractInvertedFile)

Number of indexed elements
"""
Base.length(idx::AbstractInvertedFile) = length(idx.sizes)
# SimilaritySearch.getpools(::AbstractInvertedFile, results=SimilaritySearch.GlobalKnnResult) = results
Base.show(io::IO, idx::AbstractInvertedFile) = print(io, "{$(typeof(idx)) vocsize=$(length(idx.lists)), n=$(length(idx))}")

"""
    struct InvertedFilesCaches
        Q
        P
    end

Caches used for `BinaryInvertedFile` (one per thread)

# Properties
- `Q`: posting lists involved in a query
- `P`: positions for merge algorithms
"""
struct InvertedFilesCaches
    Q::Vector{PostingList}
    P::Vector{UInt32}
end

function getcachepostinglists(pools::Vector{InvertedFilesCaches})
    Q = pools[Threads.threadid()].Q
    empty!(Q)
    Q
end

function getcachepositions(k::Integer, pools::Vector{InvertedFilesCaches})
    P = pools[Threads.threadid()].P
    resize!(P, k)
    fill!(P, 1)
    P
end

const GlobalInvertedFilesCachesPool = Vector{InvertedFilesCaches}(undef, 0)

function __init__invfile()
    n = Threads.nthreads()

    while length(GlobalInvertedFilesCachesPool) < n
        push!(GlobalInvertedFilesCachesPool, InvertedFilesCaches(Vector{PostingList}(undef, 10), Vector{UInt32}(undef, 10)))
    end

end

getpools(invfile::AbstractInvertedFile) = GlobalInvertedFilesCachesPool

"""
    sparseiterator(db, i)

Creates an iterator for indices and values of the `i`-th db's element (e.g., column).
Several specializations are provided.
"""
function sparseiterator(db::MatrixDatabase{<:SparseMatrixCSC}, i)
    r = nzrange(db.matrix, i)
    v = nonzeros(db.matrix)
    zip(r, view(v, r))
end

sparseiterator(db::MatrixDatabase{<:Matrix}, i) = enumerate(view(db.matrix, i))
sparseiterator(db::AbstractDatabase, i) = sparseiterator(db[i])

"""
    sparseiterator(obj)

`(id, weight)` iterator for `obj` for generic databases.
"""
sparseiterator(obj::KnnResult) = obj
sparseiterator(obj::AbstractVector{<:AbstractFloat}) = enumerate(obj)
sparseiterator(obj::AbstractVector{<:Integer}) = ((u, 1) for u in obj)
sparseiterator(obj::Set) = (convertpair(u) for u in obj)
sparseiterator(obj::SortedIntSet) = (convertpair(u) for u in obj)
sparseiterator(obj) = (convertpair(u) for u in obj)

"""
    convertpair(u)

Converts an element of an `sparseiterator` into an usable pair.
"""
convertpair(u::Integer) = (u, 1)
convertpair(u::Tuple) = u # assert length(u) = 2
convertpair(u::Vector) = u # assert length(u) = 2
convertpair(u::Pair) = u

"""
    Base.append!(idx, db; parallel_block=1000, pools=nothing, tol=1e-6)

Appends all `db` elements into the index `idx`. It work in parallel using all available threads.

# Arguments:
- `idx`: The inverted index
- `db`: The database of sparse objects, it can be only indices if each object is a list of integers or a set of integers (useful for `BinaryInvertedFile`),
    sparse matrices, dense matrices, among other combinations.
- `n`: The number of items to insert (defaults to all)

# Keyword arguments:
- `parallel_block`: inserts `parallel_block` elements in parallel, this argument must be larger than `Threads.nthreads()` but also not so large since the algorithm take advantage of small `parallel_block`.
- `pools`: unused argument but necessary by `searchbatch` (from `SimilaritySearch`)
- `tol`: controls what is a zero (i.e., weights < tol will be ignored).
"""
function Base.append!(idx::AbstractInvertedFile, db::AbstractDatabase, n=length(db); parallel_block=1000, pools=nothing, tol=1e-6)
    parallel_block = min(parallel_block, n)
    startID = length(idx)
    sp = 1
    resize!(idx.sizes, length(idx) + n)

    while sp < n
        ep = min(n, sp + parallel_block)
        @threads for i in sp:ep
            objID = i + startID
            nz = 0

            @inbounds for (tokenID, weight) in sparseiterator(db, i)
                weight < tol && continue
                nz += 1
                try
                    lock(idx.locks[tokenID])
                    _internal_push!(idx, tokenID, objID, weight, true)
                finally
                    unlock(idx.locks[tokenID])
                end
            end

            @inbounds idx.sizes[objID] = nz
        end

        sp = ep + 1
    end

    idx
end

"""
    push!(idx::AbstractInvertedFile, obj; pools=nothing, tol=1e-6)

Inserts a single element into the index. This operation is not thread-safe.

# Arguments
- `idx`: The inverted index
- `obj`: The object to be indexed

# Keyword arguments
- `pools`: unused argument
- `tol`: controls what is a zero (i.e., `weight < tol` will be ignored)
"""
function Base.push!(idx::AbstractInvertedFile, obj; pools=nothing, tol=1e-6)
    n = length(idx) + 1
    nz = 0

    @inbounds for (tokenID, weight) in sparseiterator(obj)
        weight < tol && continue
        nz += 1
        _internal_push!(idx, tokenID, n, weight, false)
    end

    push!(idx.sizes, nz)
    idx
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 RouteTablePropertiesFormat2 <: SwaggerModel
    routes::Any # spec type: Union{ Nothing, Vector{Route2} } # spec name: routes
    subnets::Any # spec type: Union{ Nothing, Vector{Subnet} } # spec name: subnets
    disableBgpRoutePropagation::Any # spec type: Union{ Nothing, Bool } # spec name: disableBgpRoutePropagation
    provisioningState::Any # spec type: Union{ Nothing, ProvisioningState } # spec name: provisioningState

    function RouteTablePropertiesFormat2(;routes=nothing, subnets=nothing, disableBgpRoutePropagation=nothing, provisioningState=nothing)
        o = new()
        validate_property(RouteTablePropertiesFormat2, Symbol("routes"), routes)
        setfield!(o, Symbol("routes"), routes)
        validate_property(RouteTablePropertiesFormat2, Symbol("subnets"), subnets)
        setfield!(o, Symbol("subnets"), subnets)
        validate_property(RouteTablePropertiesFormat2, Symbol("disableBgpRoutePropagation"), disableBgpRoutePropagation)
        setfield!(o, Symbol("disableBgpRoutePropagation"), disableBgpRoutePropagation)
        validate_property(RouteTablePropertiesFormat2, Symbol("provisioningState"), provisioningState)
        setfield!(o, Symbol("provisioningState"), provisioningState)
        o
    end
end # type RouteTablePropertiesFormat2

const _property_map_RouteTablePropertiesFormat2 = Dict{Symbol,Symbol}(Symbol("routes")=>Symbol("routes"), Symbol("subnets")=>Symbol("subnets"), Symbol("disableBgpRoutePropagation")=>Symbol("disableBgpRoutePropagation"), Symbol("provisioningState")=>Symbol("provisioningState"))
const _property_types_RouteTablePropertiesFormat2 = Dict{Symbol,String}(Symbol("routes")=>"Vector{Route2}", Symbol("subnets")=>"Vector{Subnet}", Symbol("disableBgpRoutePropagation")=>"Bool", Symbol("provisioningState")=>"ProvisioningState")
Base.propertynames(::Type{ RouteTablePropertiesFormat2 }) = collect(keys(_property_map_RouteTablePropertiesFormat2))
Swagger.property_type(::Type{ RouteTablePropertiesFormat2 }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_RouteTablePropertiesFormat2[name]))}
Swagger.field_name(::Type{ RouteTablePropertiesFormat2 }, property_name::Symbol) =  _property_map_RouteTablePropertiesFormat2[property_name]

function check_required(o::RouteTablePropertiesFormat2)
    true
end

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


text:
"""
	triangulate(cmd0::String="", arg1=nothing; kwargs...)

Reads one or more ASCII [or binary] files (or standard input) containing x,y[,z] and performs Delaunay
triangulation, i.e., it find how the points should be connected to give the most equilateral
triangulation possible.

Full option list at [`triangulate`]($(GMTdoc)triangulate.html)

Parameters
----------
- **C** | **slope_grid** :: [Type => Number]

    Read a slope grid (in degrees) and compute the propagated uncertainty in the
    bathymetry using the CURVE algorithm
    ($(GMTdoc)triangulate.html#a)
- **D** | **derivatives** :: [Type => Str]

    Take either the x- or y-derivatives of surface represented by the planar facets (only used when G is set).
    ($(GMTdoc)triangulate.html#d)
- **E** | **empty** :: [Type => Str | Number]

    Set the value assigned to empty nodes when G is set [NaN].
    ($(GMTdoc)triangulate.html#e)
- **G** | **grid** | **outgrid** :: [Type => Str | []]

    Use triangulation to grid the data onto an even grid (specified with R I).
    Append the name of the output grid file.
    ($(GMTdoc)triangulate.html#g)
- **I** | **inc** :: [Type => Str | Number]

    *x_inc* [and optionally *y_inc*] is the grid spacing.
    ($(GMTdoc)triangulate.html#i)
- $(GMT.opt_J)
- **M** | **network** :: [Type => Bool]

    Output triangulation network as multiple line segments separated by a segment header record.
    ($(GMTdoc)triangulate.html#m)
- **N** | **ids** :: [Type => Bool]

    Used in conjunction with G to also write the triplets of the ids of all the Delaunay vertices
    ($(GMTdoc)triangulate.html#n)
- **Q** | **voronoi** :: [Type => Str | []]

    Output the edges of the Voronoi cells instead [Default is Delaunay triangle edges]
    ($(GMTdoc)triangulate.html#q)
- $(GMT.opt_R)
- **S** | **triangles** :: [Type => Bool]

    Output triangles as polygon segments separated by a segment header record. Requires Delaunay triangulation.
    ($(GMTdoc)triangulate.html#s)
- **T** | **edges** :: [Type => Bool]

    Output edges or polygons even if gridding has been selected with the G option
    ($(GMTdoc)triangulate.html#t)
- $(GMT.opt_V)
- **Z** | **xyz** | **triplets** :: [Type => Bool]

    ($(GMTdoc)triangulate.html#z)
- $(GMT.opt_bi)
- $(GMT.opt_bo)
- $(GMT.opt_di)
- $(GMT.opt_e)
- $(GMT.opt_f)
- $(GMT.opt_h)
- $(GMT.opt_i)
- $(GMT.opt_r)
- $(GMT.opt_swap_xy)
"""
function triangulate(cmd0::String="", arg1=nothing; kwargs...)

	length(kwargs) == 0 && return monolitic("triangulate", cmd0, arg1)

	d = KW(kwargs)
	help_show_options(d)		# Check if user wants ONLY the HELP mode
	cmd, = parse_common_opts(d, "", [:R :I :V_params :bi :bo :di :e :f :h :i :r :yx])
	cmd  = parse_these_opts(cmd, d, [[:C :slope_grid], [:D :derivatives], [:E :empty], [:G :grid :outgrid], [:M :network],
                                     [:N :ids], [:S :triangles], [:T :edges], [:Z :xyz :triplets]])
	cmd = parse_Q_tri(d, [:Q :voronoi], cmd)
	if (occursin("-I", cmd) && occursin("-R", cmd) && !occursin("-G", cmd)) cmd *= " -G"  end
	if (!occursin("-G", cmd)) cmd, = parse_J(d, cmd)  end

	common_grd(d, cmd0, cmd, "triangulate ", arg1)		# Finish build cmd and run it
end

# ---------------------------------------------------------------------------------------------------
function parse_Q_tri(d::Dict, symbs::Array{Symbol}, cmd::String)
	(show_kwargs[1]) && return print_kwarg_opts(symbs, "Bool | String")	# Just print the options
	if ((val = find_in_dict(d, symbs)[1]) !== nothing)
		cmd *= " -Q";   val_ = string(val)
		(startswith(val_, "pol")) && (cmd *= "n")
    end
    return cmd
end

# ---------------------------------------------------------------------------------------------------
triangulate(arg1::Array, cmd0::String=""; kw...) = triangulate(cmd0, arg1; kw...)


text:
"""
    SyntheticDiD

Synthetic Differences-in-Differences (SDID) estimator as described in Arkhangelsky et al. (2021).
SDID is defined as

    μ̂, α̂, β̂, τ̂ = argmin ∑ₙ∑ₜ(Yᵢₜ - μ - αᵢ - βₜ - Wᵢₜτ)²ω̂ᵢλ̂ₜ

where i = 1, 2, ..., N denotes the observation units in the panel, t = 1, ..., T denotes the time
periods, Y is an (N×T) matrix of outcomes, W is an (N×T) indicator matrix of treatment status, ω is
a unit weight and λ is a time period weight.

"""
struct SyntheticDiD{T1}# <: SCM where T1
    tp::T1
    ω̂::Vector{Float64}
    λ̂::Vector{Float64}
    ŷ::Vector{Float64}
    τ̂::Vector{Float64}
    se_τ̂::Vector{Float64}
end

function SyntheticDiD(x::BalancedPanel{SingleUnitTreatment{Continuous}})
    ω̂ = zeros(length(x.is) - 1)
    λ̂ = zeros(length_T₀(x))
    ŷ = zeros(size(x.Y, 2))
    τ̂ = zeros(1)
    se_τ̂ = zeros(1)
    SyntheticDiD(x, ω̂, λ̂, ŷ, τ̂, se_τ̂)
end

function fit!(x::SyntheticDiD{BalancedPanel{SingleUnitTreatment{Continuous}}})

    (;Y, W) = x.tp
    y₁₀, y₁₁, y₀₀, y₀₁ = decompose_y(x.tp)
    T₀ = length_T₀(x.tp)
    T₁ = length_T₁(x.tp)
    Nₜᵣ = 1 # length(treated_ids(x)) - has to be one given type of BalancedPanel in this method
    Nₖₒ = size(Y, 1) - Nₜᵣ

    ### Algorithm presented in Arkhangelsky et al. (2021), v4 arXiv

    ## (1) Compute regularization parameter ζ
    # First difference of outcomes for control units before treatment time
    Δᵢₜ = Y[Not(treated_ids(x.tp)), 2:T₀] .- Y[Not(treated_ids(x.tp)), 1:T₀-1]
    Δ̄ = mean(Δᵢₜ) # Average first difference

    σ̂ = √(sum((Δᵢₜ .- Δ̄)'*(Δᵢₜ .- Δ̄)) / (Nₖₒ*(T₀ - 1)))

    ζ = (Nₜᵣ*T₁)^(1/4) * σ̂

    ## (2) Compute unit weights ω̂
    #!# Note - paper seems to generalise to multiple treated units already, here using one unit for
    # now but should be 1/Nₖₒ * sum(Y[treated_ids(x), t]) at each time period
    #!# Weights produced here are very non-sparse due to ζ*norm regularization - double check!
    function ℓᵢ(ω)
        ω₀ = ω[1]; ω = @view ω[2:end]
        (ω₀ .+ y₀₀'*ω .- y₁₀)'*(ω₀ .+ y₀₀'*ω .- y₁₀) + 1e12*(1.0 .- sum(ω))^2 + ζ^2*T₀*norm(ω)
    end

    # Bounds: ω₀ ∈ ℝ, ω ∈ Ω = {ω ∈ ℝ₊ᴺ: ∑wᵢ = 1} - essentially 0 ≤ ωᵢ ≤ 1 with sum constraint and ω₀
    # unbounded
    #!# again this will differ for multiple treated units
    lower = [-Inf; fill(0.0, Nₖₒ)]
    upper = [Inf; fill(1.0, Nₖₒ)]
    initial = [0.0; [1/Nₖₒ for _ ∈ 1:Nₖₒ]]

    res = optimize(ℓᵢ, lower, upper, initial)

    ω̂₀ = res.minimizer[1]  #!# this is probably actually unused?
    ω̂ = res.minimizer[2:end]
    x.ω̂ .= ω̂

    ## (3) Compute time weights λ̂
    # Note small regularization added as per footnote 3 of the paper
    function ℓₜ(λ)
        λ₀ = λ[1]; λ = @view λ[2:end]
        ȳ₀₁ = vec(sum(y₀₁, dims = 2) ./ size(y₀₁, 2)) # Average post-treatment outcome for control
        (λ₀ .+ y₀₀*λ .- ȳ₀₁)'*(λ₀ .+ y₀₀*λ .- ȳ₀₁) + 1e12*(1.0 .- sum(λ))^2 + ζ^2*Nₖₒ*norm(λ)
    end

    lower = [-Inf; fill(0.0, T₀)]
    upper = [Inf; fill(1.0, T₀)]
    initial = [0.0; [1/T₀ for _ ∈ 1:T₀]]

    res = optimize(ℓₜ, lower, upper, initial)

    λ̂₀ = res.minimizer[1]  #!# this is probably actually unused?
    λ̂ = res.minimizer[2:end]

    x.λ̂ .= λ̂

    ## (4) Run regression
    # Need to estimate fixed effect model with time/unit fixed effects and
    #!# To do: for covariates we want to partial them out before this regression
    # Add unit weights - 1 for treated unit, ω̂ᵢ for untreated
    ω_df = DataFrame(Symbol(x.tp.id_var) => [treated_labels(x.tp); x.tp.is[x.tp.is .!= treated_labels(x.tp)]], :ω̂ => [1.0; ω̂])
    df = leftjoin(x.tp.df, ω_df, on = x.tp.id_var)

    # Add time weights
    λ_df = DataFrame(Symbol(x.tp.t_var) => x.tp.ts, :λ̂ => [fill(1/T₁, T₁); λ̂])
    df = leftjoin(df, λ_df, on = x.tp.t_var)

    # Total weights - ωλ
    df.sdid_weight = df.ω̂ .* df.λ̂

    # Add treatment dummy based on treatment matrix
    df.W = reshape(x.tp.W', length(x.tp.ts)*length(x.tp.is))

    sdid_reg = reg(df,
        term(x.tp.outcome_var) ~ term("W") + fe(term(x.tp.id_var)) + fe(term(x.tp.t_var)),
        weights = :sdid_weight, Vcov.robust())

    # Treatment effect is estimated as the coefficient on treatment dummy
    τ̂ = coef(sdid_reg)
    x.τ̂ .= τ̂
    x.se_τ̂ .= stderror(sdid_reg)

    # Counterfactual outcome is actual outcome plus treatment effect
    ŷ₁₁ = y₁₁ .+ τ̂
    x.ŷ .= [y₁₀; ŷ₁₁]
end

#!# TO DO - show method for this type

max_stars_repo_name: nilshg/SynthControl

max_issues_repo_name: nilshg/SynthControl

max_forks_repo_name: nilshg/SynthControl


text:
# Following multiple issues over time with ordering, this attempts to have
# bunch of test cases where it's clear in what order things are tokenized / found

@testset "Ordering-1" begin
    st = raw"""
        A

        B
        """
    steps = st |> explore_md_steps
    blocks, = steps[:ocblocks]
    @test length(blocks) == 1
    @test blocks[1].name == :COMMENT
    @test isapproxstr(st |> seval, """

A

B

            """)
end

@testset "Ordering-2" begin
    st = raw"""
        A

        \begin{eqnarray}
            1 + 1 &=& 2
        \end{eqnarray}

        B
        """

    @test isapproxstr(st |> seval, raw"""

A

            \[\begin{array}{rcl}
                1 + 1 &=& 2
            \end{array}\]

B
""")
end

@testset "Ordering-3" begin
    st = raw"""
        A

        \begin{eqnarray}
            1 + 1 &=& 2
        \end{eqnarray}

        B

        C
        """
    @test isapproxstr(st |> seval, raw"""

A

            \[\begin{array}{rcl}
                1 + 1 &=& 2
            \end{array}\]

B

C
""")
end

@testset "Ordering-4" begin
    st = raw"""
        \newcommand{\eqa}[1]{\begin{eqnarray}#1\end{eqnarray}}
        A
        \eqa{
             B
        }
        C
        """

    @test isapproxstr(st |> seval, raw"""

A
            \[\begin{array}{rcl}
                B
            \end{array}\]
            C
""")
end

@testset "Ordering-5" begin
    st = raw"""
        A [❗_ongoing_ ] C
        """
    @test isapproxstr(st |> seval, raw"""

A [❗ongoing ] C

        """)
    st = raw"""
        0
        * A
            * B [❗_ongoing_ ]
        C
        """
    @test isapproxstr(st |> seval, raw"""

0

A

B [❗ongoing ]

C

            """)
end


text:
"""
    struct Preprocessor <: Any
Hold a pointer to a `clang::Preprocessor` object.
"""
struct Preprocessor
    ptr::CXPreprocessor
end

Base.unsafe_convert(::Type{CXPreprocessor}, x::Preprocessor) = x.ptr
Base.cconvert(::Type{CXPreprocessor}, x::Preprocessor) = x


text:
##########################################################################
# Copyright 2017 Samuel Ridler.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
##########################################################################

# For determining transient period duration at start and end of simulation, based on call response durations.
# Useful for deciding duration of simulation warm-up and cool-down periods.
# Loads sim from config, splits calls into batches, runs simulation for each batch,
# then saves and plots the moving average of response durations.

using JEMSS
using Statistics

# parameters
configFilename = "sim_config.xml"
outputPath = "output"
displayPlots = false

println("Loading sim")
sim = initSim(configFilename, doPrint = false)

isdir(outputPath) || mkdir(outputPath)

callArrivalTimes = [c.arrivalTime for c in sim.calls] # backup, as call arrival times will change when calls are batched

# break calls into call batches, each with callBatchNumCalls many calls, or numCallBatches many call batches
# mutates: calls
function createCallBatchesByCount!(calls::Vector{Call}, startTime::Float;
	callBatchNumCalls::Int = nullIndex, numCallBatches::Int = nullIndex)

	@assert(all(c -> c.status == callNullStatus, calls))

	numCalls = length(calls)
	@assert(sum([callBatchNumCalls, numCallBatches] .!= nullIndex) == 1)
	if callBatchNumCalls != nullIndex
		numCallBatches = Int(floor(numCalls / callBatchNumCalls))
	elseif numCallBatches != nullIndex
		callBatchNumCalls = Int(floor(numCalls / numCallBatches))
	end

	callBatches = Vector{Vector{Call}}(undef,numCallBatches)
	callBatchStartTimes = Vector{Float}(undef,numCallBatches)
	callBatchFirstCallIndex = callBatchLastCallIndex = 0 # first and last call indices for current call batch
	for callBatchIndex = 1:numCallBatches
		# get call batch start time
		if callBatchIndex == 1
			callBatchStartTimes[1] = startTime
		else
			callBatchStartTimes[callBatchIndex] = calls[callBatchLastCallIndex].arrivalTime
		end
		# find first and last calls for call batch callBatchIndex
		callBatchFirstCallIndex = callBatchLastCallIndex + 1
		callBatchLastCallIndex = callBatchFirstCallIndex + callBatchNumCalls - 1
		callBatches[callBatchIndex] = calls[callBatchFirstCallIndex:callBatchLastCallIndex]
	end

	# renumber calls from 1 to n within each call batch, set each call batch to start at startTime
	for (callBatchIndex, callBatch) in enumerate(callBatches)
		for (callIndex, call) in enumerate(callBatch)
			call.index = callIndex
			call.arrivalTime -= (callBatchStartTimes[callBatchIndex] - startTime)
		end
	end

	return callBatches
end

# run simulation for each callBatch in callBatches
# mutates: sim, callBatches
function simulateCallBatches!(sim::Simulation, callBatches::Vector{Vector{Call}})
	# run simulation for each call batch
	numCallBatches = length(callBatches)
	for (callBatchIndex, callBatch) in enumerate(callBatches)
		print("\rSimulating call batch ", callBatchIndex, " of ", numCallBatches)
		setSimCalls!(sim, callBatch)
		simulate!(sim)
		sim.calls = Call[] # remove call batch
	end
	println()
end

function calcAvgResponseDurations(callBatches::Vector{Vector{Call}})
	# for callBatches that have already been used in simulation, calculate avgResponseDurations,
	# where avgResponseDurations[i] is the average response duration of call i from all batches
	@assert(all(callBatch -> all(call -> call.responseDuration != nullTime, callBatch), callBatches))
	avgResponseDurations = [mean([callBatches[i][j].responseDuration for i = 1:length(callBatches)]) for j = 1:length(callBatches[1])] * (24*60)
	return avgResponseDurations
end

function calcMovingAvgResponseDurations(avgResponseDurations::Vector{Float}, movingWindowSize::Int)
	callBatchNumCalls = length(avgResponseDurations)

	# # for constant window size:
	# numWindows = max(0, callBatchNumCalls - movingWindowSize + 1)
	# movingAvgResponseDurations = [mean(avgResponseDurations[i + (0:movingWindowSize-1)]) for i = 1:numWindows]

	# # for window size = min(window size, number of calls so far):
	# numWindows = callBatchNumCalls
	# movingAvgResponseDurations = [mean(avgResponseDurations[max(1,i-(movingWindowSize-1)):i]) for i = 1:numWindows]

	# for window size = 1 + 2 * min(movingWindowSize, max allowable window size):
	# (each average value is at centre of window)
	numWindows = callBatchNumCalls
	movingAvgResponseDurations = Vector{Float}(undef,numWindows)
	for windowIndex = 1:numWindows
		currentWindowSize = min(movingWindowSize, windowIndex - 1, numWindows - windowIndex)
		movingAvgResponseDurations[windowIndex] = mean(avgResponseDurations[windowIndex .+ (-currentWindowSize : currentWindowSize)])
	end

	# # for window size = 1 + 2 * min(movingWindowSize, max allowable window size):
	# # (each average value is at centre of window, last n (=movingWindowSize) average values are not calculated)
	# # this follows Welch's procedure in "On the problem of the initial transient in steady-state simulation"
	# numWindows = max(0, callBatchNumCalls - movingWindowSize)
	# movingAvgResponseDurations = Vector{Float}(undef,numWindows)
	# for windowIndex = 1:numWindows
		# currentWindowSize = min(movingWindowSize, max(0, windowIndex-1))
		# movingAvgResponseDurations[windowIndex] = mean(avgResponseDurations[windowIndex .+ (-currentWindowSize : currentWindowSize)])
	# end

	return movingAvgResponseDurations
end

# For calls in sim.calls, split into n batches (n = callSetsNumCallBatches[i], i=1,...), then simulate.
# The response duration of the jth call is averaged over all batches.
function getCallSetsAvgResponseDurations(sim, callSetsNumCallBatches::Vector{Int})::Vector{Vector{Float}}
	numCallSets = length(callSetsNumCallBatches)
	callSetsAvgResponseDurations = []
	calls = sim.calls # shorthand
	for (callSetIndex, numCallBatches) in enumerate(callSetsNumCallBatches)
		println("Call set ", callSetIndex, " of ", numCallSets)

		# create batches (for current set), simulate for each, calculate average response durations
		reset!(sim)
		callBatches = createCallBatchesByCount!(calls, sim.startTime; numCallBatches = numCallBatches)
		simulateCallBatches!(sim, callBatches)
		push!(callSetsAvgResponseDurations, calcAvgResponseDurations(callBatches))

		# undo changes from createCallBatchesByCount!(), and reset calls
		for (i,c) in enumerate(calls)
			c.index = i
			c.arrivalTime = callArrivalTimes[i]
		end
		setSimCalls!(sim, calls)
	end

	return callSetsAvgResponseDurations
end

function writeCallSetsAvgResponseDurationsFile(filename::String, callSetsAvgResponseDurations::Vector{Vector{Float}}, callSetsNumCallBatches::Vector{Int}, movingWindowSizes::Vector{Int})
	tables = Table[]
	for (i, numCallBatches) in enumerate(callSetsNumCallBatches)
		avgResponseDurations = callSetsAvgResponseDurations[i]
		m = length(avgResponseDurations)
		header = ["callIndex"]
		data = Array{Any,2}(undef, m, length(movingWindowSizes) + 1) # to be filled
		data[:] .= "" # empty values
		data[:,1] = 1:m
		for (j, movingWindowSize) in enumerate(movingWindowSizes)
			push!(header, "windowSize_$movingWindowSize")
			movingAvgResponseDurations = calcMovingAvgResponseDurations(avgResponseDurations, movingWindowSize)
			if movingAvgResponseDurations != []
				data[1:length(movingAvgResponseDurations), j+1] = movingAvgResponseDurations
			end
		end
		push!(tables, Table("numCallBatches_$numCallBatches", header, data))
	end
	writeTablesToFile(filename, tables)
end

# Given the response duration of the jth call averaged over all batches, moving averages are calculated using
# each value in movingWindowSizes. One plot is produced for each value in callSetsNumCallBatches.
function plotMovingAvgResponseDurations(callSetsAvgResponseDurations::Vector{Vector{Float}}, callSetsNumCallBatches::Vector{Int}, movingWindowSizes::Vector{Int})
	numCallSets = length(callSetsNumCallBatches)
	@assert(numCallSets == length(callSetsAvgResponseDurations))
	plots = []
	for i = 1:numCallSets
		# create plot of moving average response durations with different moving window sizes
		push!(plots, Plots.plot(title = string("Number of call batches: ", callSetsNumCallBatches[i]),
			xaxis = "Call index", yaxis = "Mean response duration (minutes)"))
		for movingWindowSize in movingWindowSizes
			movingAvgResponseDurations = calcMovingAvgResponseDurations(callSetsAvgResponseDurations[i], movingWindowSize)
			if movingAvgResponseDurations != []
				Plots.plot!(plots[i], movingAvgResponseDurations, label = string(movingWindowSize));
			end
		end
	end
	return plots
end

println("Simulating with different call batches")
callSetsNumCallBatches = [5,10] # for runtests.jl
# callSetsNumCallBatches = [100:100:1000;] # for actual use
movingWindowSizes = [0, 1, 2, 5, 10, 20, 50, 100, 200, 500, 1000];
callSetsAvgResponseDurations = getCallSetsAvgResponseDurations(sim, callSetsNumCallBatches)

println("Saving results to file")
filename = joinpath(outputPath, "transient.csv")
writeCallSetsAvgResponseDurationsFile(filename, callSetsAvgResponseDurations, callSetsNumCallBatches, movingWindowSizes)

if displayPlots
	println("Plotting")
	import Plots
	Plots.plotly()
	plots = plotMovingAvgResponseDurations(callSetsAvgResponseDurations, callSetsNumCallBatches, movingWindowSizes)
	for plot in plots display(plot) end
end


text:
println("<HEARTBEAT>")
using Test, Mill, Flux
using Random

function ngradient(f, xs::AbstractArray...)
  grads = zero.(xs)
  for (x, Δ) in zip(xs, grads), i in 1:length(x)
    δ = sqrt(eps())
    tmp = x[i]
    x[i] = tmp - δ/2
    y1 = f(xs...)
    x[i] = tmp + δ/2
    y2 = f(xs...)
    x[i] = tmp
    Δ[i] = (y2-y1)/δ
  end
  return grads
end

gradcheck(f, xs...) =
  all(isapprox.(ngradient(f, xs...),
                gradient(f, xs...), rtol = 1e-5, atol = 1e-5))

gradtest(f, xs::AbstractArray...) = gradcheck((xs...) -> sum(sin.(f(xs...))), xs...)
gradtest(f, dims...) = gradtest(f, rand.(Float64, dims)...)

function mngradient(f, xs::AbstractArray...)
    grads = zero.(xs)
    for (x, Δ) in zip(xs, grads), i in eachindex(x)
        δ = sqrt(eps())
        tmp = x[i]
        x[i] = tmp - δ/2
        y1 = f(xs...)
        x[i] = tmp + δ/2
        y2 = f(xs...)
        x[i] = tmp
        Δ[i] = (y2-y1)/δ
    end
    return grads
end

function mgradcheck(f, xs...)
    correct = true
    for (ng, ag) in zip(mngradient(f, xs...), gradient(f, xs...))
        if isnothing(ag)
            if any(ng .!= 0)
                correct = false
                @show ng
            end
        elseif !all(isapprox.(ng, ag, rtol = 1e-5, atol = 1e-5))
            correct = false
            grad_dif = [abs.(x) for x in (ng .- ag)]
            @show ng
            @show ag
            @show grad_dif
        end
    end
    correct
end

mgradtest(f, xs::AbstractArray...) = mgradcheck((xs...) -> sum(sin.(f(xs...))), xs...)

Random.seed!(24)

const BAGS = [
      length2bags([1, 2, 3]),
      length2bags([4, 1, 1]),
      length2bags([2, 2, 2])
]

const BAGS2 = [
        length2bags(ones(Int, 10)),
        length2bags(2 .* ones(Int, 5)),
        length2bags([5, 5]),
        length2bags([10]),
        length2bags([3, 4, 3]),
        AlignedBags([1:3, 0:-1, 0:-1, 4:7, 0:-1, 8:10]),
        AlignedBags([0:-1, 1:5, 0:-1, 0:-1, 0:-1, 6:10]),
        ScatteredBags([collect(1:3), collect(7:10), collect(4:6)]),
        ScatteredBags([collect(7:10), [], collect(1:3), [], collect(4:6), []]),
        ScatteredBags([[], collect(1:10), []]),
]

const BAGS3 = [
         (AlignedBags([1:2, 3:4, 0:-1]), ScatteredBags([[2,3,4], [1], []]), AlignedBags([1:4, 0:-1, 5:8, 0:-1])),
         (AlignedBags([0:-1, 1:2, 3:4]), ScatteredBags([[1], [2], [3, 4]]), AlignedBags([0:-1, 1:7, 0:-1, 8:8])),
         (AlignedBags([0:-1, 0:-1, 1:2, 3:4]), ScatteredBags([[2,4], [], [3, 1], []]), AlignedBags([1:1, 2:2, 0:-1, 3:8])),
         (AlignedBags([0:-1, 1:2, 3:4, 0:-1]), ScatteredBags([[], [1,3], [2,4], []]), AlignedBags([0:-1, 1:2, 3:6, 7:8]))
        ]

println("<HEARTBEAT>")

@testset "Data nodes" begin
    include("datanode.jl")
end
@testset "Model nodes" begin
    include("modelnode.jl")
end
@testset "Missing" begin
    include("missing.jl")
end
println("<HEARTBEAT>")
@testset "Aggregation" begin
    include("aggregation.jl")
end
println("<HEARTBEAT>")
@testset "Gradtests" begin
    include("gradtests.jl")
end
println("<HEARTBEAT>")
@testset "Conv" begin
    include("conv.jl")
end
@testset "Bags" begin
    include("bags.jl")
end
@testset "NGrams" begin
    include("ngrams.jl")
end
@testset "Activations" begin
    include("activations.jl")
end
@testset "Hierarchical Utils" begin
    include("hierarchical_utils.jl")
end
@testset "Replace in" begin
    include("replacein.jl")
end
@testset "Partial Eval" begin
    include("partialeval.jl")
end
@testset "Lazy Node" begin
    include("lazynode.jl")
end


text:
mutable struct Variable{T,U <: Union{Vector{T},Optional{T},T}}
    value::U
    default_value::U
end

function get_value(v::Variable{T,Vector{T}}) where T
    !isempty(v.value) && return v.value
    return v.default_value
end

function get_value(v::Variable)
    v.value !== nothing && return v.value
    return v.default_value
end

set_value(v::Variable{T,Vector{T}}, arg) where T = push!(v.value, arg)
set_value(v::Variable, arg) = (v.value = arg)

function empty_var(arg::Argument{T,<:Any,V}) where {T,V}
    dest_is_vector = V
    if dest_is_vector
        default = arg.default !== nothing ? arg.default : T[]
        value = arg.action in [:append, :append_const] ? copy(default) : T[]
        return Variable{T,Vector{T}}(value, default)
    end

    if arg.action === :count
        default = arg.default !== nothing ? arg.default : 0
        return Variable{T,T}(arg.default, arg.default)
    end

    return Variable{T,Optional{T}}(nothing, arg.default)
end

function promote_var(arg1::Variable{S,V1}, arg2::Variable{T,V2}) where {S,T,V1,V2}
    if ((V1 <: AbstractVector && !(V2 <: AbstractVector)) ||
        (V2 <: AbstractVector && !(V1 <: AbstractVector)))

        throw(ArgumentError("Multiple arguments target the same destination variable,\n" *
        "but one accepts multiple values, and the other only accepts single values"))
    end

    if (arg1.default_value !== nothing && arg2.default_value !== nothing &&
        arg1.default_value != arg2.default_value)
        throw(ArgumentError("Multiple arguments target the same destination variable,\n" *
        "but have different default values(`$(arg1.default_value)` != `$(arg2.default_value)`"))
    end

    default_value = something(arg1.default_value, Some(arg2.default_value))

    U = promote_type(S, T)
    dest_is_vector = V1 <: AbstractVector

    if dest_is_vector
        return Variable{U,Vector{U}}(U[], default_value)
    end

    return Variable{U,Optional{U}}(nothing, default_value)
end


text:
module NotebookToLaTeX

using Requires
using ReadableRegex
# using Plots
# using Makie
# using CairoMakie
using Base64
using JSON
using Librsvg_jll # For converting svg to pdf

export notebooktolatex, jupytertolatex

include("templates.jl")
include("auxiliarytex.jl")
include("markdowntolatex.jl")
include("helperfunctions.jl")
# include("plutoconverter.jl")
include("jupyterconverter.jl")

export nestedget

"""
    Runner is a module for controling the scope
    when running the notebook files, avoinding that any
    command interferes with the "outside" script.
"""
module Runner
end

function __init__()

    @require Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" begin
    @require Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" begin
        @require CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0" include("plutoconverter.jl")
    end
    end
end

"""
    notebooktolatex(notebook, targetdir="./build_latex"; template=:book, fontpath=nothing)
Takes a notebook file, converts it to LaTeX and creates a file structure
with figures, fonts and listing files.
* `targetdir` is the target directory where the LaTeX project will be created.
If the directory does not exist, it is created.
* `template` - The template for the LaTeX file. It's based on LaTeX templates.
Current supported templates are `:book`, `:mathbook`.
* `fontpath` - The output LaTeX files uses JuliaMono fonts in order to support the
unicodes that are also supported in Julia. If the user already has JuliaMono installed,
he can provide the path to where the `.ttf` files are stored. If `nothing` is passed,
then the font files will be downloaded and saved in the `./fonts/` folder.
"""
function notebooktolatex(notebook::String, targetdir="./build_latex"; template=:book, fontpath=nothing)
    if endswith((notebook),".jl")
        plutotolatex(notebook, targetdir, template=template, fontpath=fontpath)
    elseif endswith((notebook),".ipynb")
        jupytertolatex(notebook, targetdir, template=template, fontpath=fontpath)
    else
        throw(ArgumentError(notebook, "extension must be either .jl or .ipynb"))
    end
end

function dispatch_output(command_eval, notebookname, path, figureindex)
   return command_eval
end

end


text:
using
  Documenter,
  SeawaterPolynomials

format = Documenter.HTML(
    collapselevel = 2,
       prettyurls = get(ENV, "CI", nothing) == "true",
        canonical = "https://clima.github.io/SeawaterPolynomials/stable/",
       mathengine = MathJax3()
)

makedocs(
    modules = [SeawaterPolynomials],
    doctest = false,
      clean = true,
   checkdocs = :all,
     format = format,
    authors = "Gregory L. Wagner and Ali Ramadhan",
   sitename = "SeawaterPolynomials.jl",

      pages = Any[
              "Home" => "index.md",
              "API" => "API.md",
                 ]
)

deploydocs(
        repo = "github.com/CliMA/SeawaterPolynomials.jl.git",
    versions = ["stable" => "v^", "v#.#", "dev" => "dev"],
push_preview = true
)


text:
using Documenter, SFGTools

makedocs()


text:
using DashBootstrapComponents, DashHtmlComponents

tabs = html_div([
    dbc_tabs([
        dbc_tab(label="Tab 1", tab_style=Dict("marginLeft" => "auto")),
        dbc_tab(label="Tab 2", label_style=Dict("color" => "#00AEF9")),
    ]),
    html_br(),
    dbc_tabs([
        dbc_tab(label="Tab 1", tabClassName="ml-auto"),
        dbc_tab(label="Tab 2", labelClassName="text-success"),
    ]),
]);

max_stars_repo_name: harisbal/dash-bootstrap-components

max_issues_repo_name: harisbal/dash-bootstrap-components

max_forks_repo_name: harisbal/dash-bootstrap-components


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder

name = "texinfo"
version = v"6.5.0"

# Collection of sources required to build texinfo
sources = [
    "https://ftp.gnu.org/gnu/texinfo/texinfo-6.5.tar.xz" =>
    "77774b3f4a06c20705cc2ef1c804864422e3cf95235e965b1f00a46df7da5f62",

]

# Bash recipe for building across all platforms
script = raw"""
cd $WORKSPACE/srcdir
cd texinfo-6.5/
./configure --prefix=$prefix --host=$target
make
make install
exit

"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = [
    Linux(:x86_64, libc=:glibc)
]

# The products that we will ensure are always built
products(prefix) = [
    ExecutableProduct(prefix, "makeinfo", Symbol("makeinfo")),
    ExecutableProduct(prefix, "pdftexi2dvi", Symbol("pdftexi2dvi")),
    ExecutableProduct(prefix, "pod2texi", Symbol("pod2texi")),
    ExecutableProduct(prefix, "texi2any", Symbol("texi2any")),
    ExecutableProduct(prefix, "texi2dvi", Symbol("texi2dvi")),
    ExecutableProduct(prefix, "texi2pdf", Symbol("texi2pdf")),
    ExecutableProduct(prefix, "texindex", Symbol("texindex"))
]

# Dependencies that must be installed before this package can be built
dependencies = [

]

# Build the tarballs, and possibly a `build.jl` as well.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies)

max_stars_repo_name: yukota/TexinfoBuilder

max_issues_repo_name: yukota/TexinfoBuilder

max_forks_repo_name: yukota/TexinfoBuilder


text:
# auxiliary functions:
# calculate plotting coordinates of nodes & edges for PhyloNetworks objects

"""
    getEdgeNodeCoordinates(net, useEdgeLength)

Calculate coordinates for plotting later with Gadfly or RCall.

Actually modifies some (minor) attributes of the network,
as it calls `directEdges!`, `preorder!` and `cladewiseorder!`.
"""
function getEdgeNodeCoordinates(net::HybridNetwork, useEdgeLength::Bool)
    try
        directEdges!(net)   # to update isChild1
    catch e
        if isa(e, PhyloNetworks.RootMismatch)
            e = PhyloNetworks.RootMismatch( e.msg * "\nPlease change the root, perhaps using rootatnode! or rootatedge!")
        end
        rethrow(e)
    end
    preorder!(net)       # to update net.nodes_changed: true pre-ordering
    cladewiseorder!(net) # to update cladewiseorder_nodeIndex: cladewise on major tree

    # determine y for each node = y of its parent edge: post-order traversal
    # also [yB,yE] for each internal node: range of y's of all children nodes
    ymin = 1.0; ymax = Float64(net.numTaxa);
    node_y  = zeros(Float64, net.numNodes) # order: in net.nodes, *!not in nodes_changed!*
    node_yB = zeros(Float64,net.numNodes) # min (B=begin) and max (E=end)
    node_yE = zeros(Float64,net.numNodes) #   of at children's nodes
    # set node_y of leaves: follow cladewise order
    nexty = ymax # first tips at the top, last at bottom
    for i=length(net.node):-1:1
        ni = net.cladewiseorder_nodeIndex[i]
        if net.node[ni].leaf
            node_y[ni] = nexty
            nexty -= 1.0
        end
    end
    # set node_y of internal nodes: follow post-order
    for i=length(net.node):-1:1
        nn = net.nodes_changed[i]
        !nn.leaf || continue # previous loop took care of leaves
        ni = findfirst(x -> x===nn, net.node)
        node_yB[ni]=ymax; node_yE[ni]=ymin;
        minor_yB  = ymax; minor_yE  = ymin;
        nomajorchild=true
        for e in nn.edge
            if nn == PhyloNetworks.getParent(e) # if e = child of node
                if e.isMajor || nomajorchild
                    cc = PhyloNetworks.getChild(e)
                    yy = node_y[findfirst(x -> x===cc, net.node)]
                    yy!==nothing || error("oops, child $(cc.number) has not been visited before node $(nn.number).")
                end
                if e.isMajor
                    nomajorchild = false # we found a child edge that is a major edge
                    node_yB[ni] = min(node_yB[ni], yy)
                    node_yE[ni] = max(node_yE[ni], yy)
                elseif nomajorchild # e is minor edge, and no major found so far
                    minor_yB = min(minor_yB, yy)
                    minor_yE = max(minor_yE, yy)
                end
            end
        end
        if nomajorchild # if all children edges are minor hybrid edges. If so: not level 1, not tree-child
            if minor_yB == minor_yE # one single child. jitter by 0.1 to make the plot readable
                minor_yB += (minor_yB < (ymax+ymin)/2 ? 0.1 : -0.1)
                minor_yE = minor_yB
            end
            node_yB[ni] = minor_yB
            node_yE[ni] = minor_yE
        end
        node_y[ni] = (node_yB[ni]+node_yE[ni])/2
    end

    # setting branch lengths for plotting
    elenCalculate = !useEdgeLength
    if (useEdgeLength)
        allBLmissing = true; nonBLmissing = true;
        for e in net.edge
            if (nonBLmissing && e.length==-1.0) nonBLmissing=false; end
            if (allBLmissing && e.length!=-1.0) allBLmissing=false; end
        end
        if (allBLmissing)
            println("All edge lengths are missing, won't be used for plotting.")
            elenCalculate = true
        end
        if (!nonBLmissing && !allBLmissing) # not all, but some are missing
            @warn "At least one non-missing edge length: plotting any missing length as 1.0"
        end
    end
    elen = Float64[] # edge lengths to be used for plotting. same order as net.edge.
    if (elenCalculate)
        # setting elen such that the age of each node = 1 + age of oldest child
        # (including minor hybrid edges): need true post-ordering.
        # calculating node ages first, elen will be calculated later.
        elen     = zeros(Float64,net.numEdges)
        node_age = zeros(Float64,net.numNodes)
        for i=length(net.node):-1:1 # post-order traversal
            if (net.nodes_changed[i].leaf) continue; end
            ni = findfirst(x -> x===net.nodes_changed[i], net.node)
            for e in net.nodes_changed[i].edge # loop over children only
                if net.nodes_changed[i] == (e.isChild1 ? e.node[2] : e.node[1])
                    node_age[ni] = max(node_age[ni], 1 +
                     node_age[findfirst(x -> x=== PhyloNetworks.getChild(e), net.node)])
                end
            end
        end
    else
        for e in net.edge
            push!(elen, (e.length==-1.0 ? 1.0 : e.length))
        end
    end

    # determine xB,xE for each edge: pre-order traversal, uses branch lengths
    # then x and yB,yE for each node: x=xE of parent edge
    xmin = 1.0; xmax=xmin
    node_x  = zeros(Float64,net.numNodes) # order: in net.nodes, *!not in nodes_changed!*
    edge_xB = zeros(Float64,net.numEdges) # min (B=begin) and max (E=end)
    edge_xE = zeros(Float64,net.numEdges) # xE-xB = edge length
    edge_yE = zeros(Float64,net.numEdges) # yE of edge = y of child node
    node_x[net.root] = xmin # root node: x=xmin=0
    for i=2:length(net.node)              # true pre-order, skipping the root (i=1)
        ni = findfirst(x -> x===net.nodes_changed[i], net.node)
        ei = nothing # index of major parent edge of current node
        for e in net.nodes_changed[i].edge
            if (e.isMajor && net.nodes_changed[i] == e.node[e.isChild1 ? 1 : 2]) # major parent edge
                ei = findfirst(x -> x===e, net.edge)
                break
            end
        end
        ei !== nothing || error("oops, could not find major parent edge of node number $ni.")
        edge_yE[ei] = node_y[ni]
        pni = findfirst(x -> x===PhyloNetworks.getParent(net.edge[ei]), net.node) # parent node index
        edge_xB[ei] = node_x[pni]
        if elenCalculate
            elen[ei] = node_age[pni] - node_age[ni]
        end
        edge_xE[ei] = edge_xB[ei] + elen[ei]
        node_x[ni] = edge_xE[ei]
        xmax = max(xmax, edge_xE[ei])
    end
    edge_yB = copy(edge_yE) # true for tree and major edges
    for i=1:net.numEdges
        if (!net.edge[i].isMajor) # minor hybrid edges
            cni = findfirst(x -> x===PhyloNetworks.getChild( net.edge[i]), net.node)
            pni = findfirst(x -> x===PhyloNetworks.getParent(net.edge[i]), net.node)
            # indices of child and parent nodes
            edge_xB[i] = node_x[pni]
            edge_xE[i] = node_x[cni]
            edge_yB[i] = node_y[pni]
            edge_yE[i] = node_y[cni]
            #@show i; @show net.edge[i]; @show pni; @show net.node[pni]; @show cni; @show net.node[cni]
        end
    end

    #@show node_x;  @show node_yB; @show node_y;  @show node_yE
    #@show edge_xB; @show edge_xE; @show edge_yB; @show edge_yE
    return edge_xB, edge_xE, edge_yB, edge_yE,
           node_x, node_y, node_yB, node_yE,
           xmin, xmax, ymin, ymax
end

"""
    checkNodeDataFrame(net, nodeLabel)

Check data frame for node annotations:
- check that the data has at least 2 columns (if it has any)
- check that the first column has integers (to serve as node numbers)
- remove rows with no node numbers
- warning if some node numbers in the data are not in the network.
"""
function checkNodeDataFrame(net::HybridNetwork, nodeLabel::DataFrame)
    labelnodes = size(nodeLabel,1)>0
    if (labelnodes && (size(nodeLabel,2)<2 || !(nonmissingtype(eltypes(nodeLabel)[1]) <: Integer)))
        @warn "nodeLabel should have 2+ columns, the first one giving the node numbers (Integer)"
        labelnodes = false
    end
    if labelnodes # remove rows with no node number, check if at least one row remains
        filter!(row->!ismissing(row[1]), nodeLabel)
        labelnodes = size(nodeLabel,1)>0
    end
    if labelnodes
      tmp = setdiff(nodeLabel[!,1], [n.number for n in net.node])
      if length(tmp)>0
        msg = "Some node numbers in the nodeLabel data frame are not found in the network:\n"
        for a in tmp msg *= string(" ",a); end
        @warn msg
      end
    end
    return(labelnodes, nodeLabel)
end

"""
    prepareNodeDataFrame(net, nodeLabel::DataFrame,
        showNodeNumber::Bool, showIntNodeLabel::Bool, labelnodes::Bool,
        node_x, node_y)

Make data frame for node annotation. `node_*` should be Float64 vectors.
`nodeLabel` should have columns as required by [`checkNodeDataFrame`](@ref)

Columns of output data frame:
- x, y: coordinates on the plots (from `node_*`)
- name: node name
- num: node number
- lab: node label
- lea: is leaf?
"""
function prepareNodeDataFrame(net::HybridNetwork, nodeLabel::DataFrame,
        showNodeNumber::Bool, showIntNodeLabel::Bool, labelnodes::Bool,
        node_x::Array{Float64,1}, node_y::Array{Float64,1})
    nrows = (showNodeNumber || showIntNodeLabel || labelnodes ? net.numNodes : net.numTaxa)
    ndf = DataFrame([String,String,String,Bool,Float64,Float64], # column types, column names, nrows
           [Symbol("name"),Symbol("num"),Symbol("lab"),Symbol("lea"),Symbol("x"),Symbol("y")], nrows)
    j=1
    for i=1:net.numNodes
    if net.node[i].leaf  || showNodeNumber || showIntNodeLabel || labelnodes
        ndf[j,:name] = net.node[i].name
        ndf[j,:num] = string(net.node[i].number)
        if labelnodes
          jn = findfirst(isequal(net.node[i].number), nodeLabel[!,1])
          ndf[j,:lab] = (jn===nothing || ismissing(nodeLabel[jn,2]) ? "" :  # node label not in table or missing
            (nonmissingtype(eltypes(nodeLabel)[2]) <: AbstractFloat ?
              @sprintf("%0.3g",nodeLabel[jn,2]) : string(nodeLabel[jn,2])))
        end
        ndf[j,:lea] = net.node[i].leaf # use this later to remove H? labels
        ndf[j,:y] = node_y[i]
        ndf[j,:x] = node_x[i]
        j += 1
    end
    end
    # @show ndf
    return(ndf)
end

"""
    prepareEdgeDataFrame(net, edgeLabel::DataFrame, mainTree::Bool,
        edge_xB, edge_xE, edge_yB, edge_yE)

Check data frame for edge annotation.
`edge_*`: Float64 vectors giving the coordinates for the beginning and end of edges.
Return data frame with columns
- x, y: coordinates on the plots
- len: node name
- gam: gamma (inheritance value)
- num: node number
- lab: node label
- hyb: is hybrid?
- min: is minor?
"""
function prepareEdgeDataFrame(net::HybridNetwork, edgeLabel::DataFrame, mainTree::Bool,
        edge_xB::Array{Float64,1}, edge_xE::Array{Float64,1},
        edge_yB::Array{Float64,1}, edge_yE::Array{Float64,1})
    nrows = net.numEdges - (mainTree ? net.numHybrids : 0)
    edf = DataFrame([String,String,String,String,Bool,Bool,Float64,Float64],
                  [Symbol("len"),Symbol("gam"),Symbol("num"),Symbol("lab"),
                   Symbol("hyb"),Symbol("min"),Symbol("x"),Symbol("y")], nrows)
    labeledges = size(edgeLabel,1)>0
    if (labeledges && (size(edgeLabel,2)<2 || !(nonmissingtype(eltypes(edgeLabel)[1]) <: Integer)))
        @warn "edgeLabel should have 2+ columns, the first one giving the edge numbers (Integer)"
        labeledges = false
    end
    if labeledges # remove rows with no edge number and check if at least one remains
        filter!(row->!ismissing(row[1]), edgeLabel)
        labeledges = size(edgeLabel,1)>0
    end
    if labeledges
      tmp = setdiff(edgeLabel[!,1], [e.number for e in net.edge])
      if length(tmp)>0
        msg = "Some edge numbers in the edgeLabel data frame are not found in the network:\n"
        for a in tmp msg *= string(" ",a); end
        @warn msg
      end
    end
    j=1
    for i = 1:length(net.edge)
        if (!mainTree || !net.edge[i].hybrid || net.edge[i].isMajor)
            edf[j,:len] = (net.edge[i].length==-1.0 ? "" : @sprintf("%0.3g",net.edge[i].length))
            # @sprintf("%c=%0.3g",'γ',net.edge[i].length)
            edf[j,:gam] = (net.edge[i].gamma==-1.0  ? "" : @sprintf("%0.3g",net.edge[i].gamma))
            edf[j,:num] = string(net.edge[i].number)
            if labeledges
              je = findfirst(isequal(net.edge[i].number), edgeLabel[!,1])
              edf[j,:lab] = (je===nothing || ismissing(edgeLabel[je,2]) ? "" :  # edge label not found in table
                (nonmissingtype(eltypes(edgeLabel)[2]) <: AbstractFloat ?
                  @sprintf("%0.3g",edgeLabel[je,2]) : string(edgeLabel[je,2])))
            end
            edf[j,:hyb] = net.edge[i].hybrid
            edf[j,:min] = !net.edge[i].isMajor
            edf[j,:y] = (edge_yB[i] + edge_yE[i])/2
            edf[j,:x] = (edge_xB[i] + edge_xE[i])/2
            j += 1
        end
    end
    # @show edf
    return labeledges, edf
end


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

@doc raw"""UserInfo holds the information about the user needed to implement the user.Info interface.

    IoK8sApiAuthenticationV1UserInfo(;
        extra=nothing,
        groups=nothing,
        uid=nothing,
        username=nothing,
    )

    - extra::Dict{String, Vector{String}} : Any additional information provided by the authenticator.
    - groups::Vector{String} : The names of groups this user is a part of.
    - uid::String : A unique value that identifies this user across time. If this user is deleted and another user by the same name is added, they will have different UIDs.
    - username::String : The name that uniquely identifies this user among all active users.
"""
mutable struct IoK8sApiAuthenticationV1UserInfo <: SwaggerModel
    extra::Any # spec type: Union{ Nothing, Dict{String, Vector{String}} } # spec name: extra
    groups::Any # spec type: Union{ Nothing, Vector{String} } # spec name: groups
    uid::Any # spec type: Union{ Nothing, String } # spec name: uid
    username::Any # spec type: Union{ Nothing, String } # spec name: username

    function IoK8sApiAuthenticationV1UserInfo(;extra=nothing, groups=nothing, uid=nothing, username=nothing)
        o = new()
        validate_property(IoK8sApiAuthenticationV1UserInfo, Symbol("extra"), extra)
        setfield!(o, Symbol("extra"), extra)
        validate_property(IoK8sApiAuthenticationV1UserInfo, Symbol("groups"), groups)
        setfield!(o, Symbol("groups"), groups)
        validate_property(IoK8sApiAuthenticationV1UserInfo, Symbol("uid"), uid)
        setfield!(o, Symbol("uid"), uid)
        validate_property(IoK8sApiAuthenticationV1UserInfo, Symbol("username"), username)
        setfield!(o, Symbol("username"), username)
        o
    end
end # type IoK8sApiAuthenticationV1UserInfo

const _property_map_IoK8sApiAuthenticationV1UserInfo = Dict{Symbol,Symbol}(Symbol("extra")=>Symbol("extra"), Symbol("groups")=>Symbol("groups"), Symbol("uid")=>Symbol("uid"), Symbol("username")=>Symbol("username"))
const _property_types_IoK8sApiAuthenticationV1UserInfo = Dict{Symbol,String}(Symbol("extra")=>"Dict{String, Vector{String}}", Symbol("groups")=>"Vector{String}", Symbol("uid")=>"String", Symbol("username")=>"String")
Base.propertynames(::Type{ IoK8sApiAuthenticationV1UserInfo }) = collect(keys(_property_map_IoK8sApiAuthenticationV1UserInfo))
Swagger.property_type(::Type{ IoK8sApiAuthenticationV1UserInfo }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiAuthenticationV1UserInfo[name]))}
Swagger.field_name(::Type{ IoK8sApiAuthenticationV1UserInfo }, property_name::Symbol) =  _property_map_IoK8sApiAuthenticationV1UserInfo[property_name]

function check_required(o::IoK8sApiAuthenticationV1UserInfo)
    true
end

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


text:
#=
    A very simple example of a semidefinite program
    in standard form. The problem is given by
    min tr(CX)
    st  tr(AX) = b
        X ⪰ 0
=#

using MySDPSolver
using JuMP
using LinearAlgebra

# Data
A = [1.0 0.0
     0.0 0.0]
b = 1.0
C = [2.0 1.0
     1.0 2.0]
# Set the model
model = Model(MySDPSolver.Optimizer)
@variable(model, X[1:2,1:2], PSD)
@constraint(model, tr(A*X) == b)
@objective(model, Min, tr(C*X))
optimize!(model)

max_stars_repo_name: TomaszWasTaken/MySDPSolver

max_issues_repo_name: TomaszWasTaken/MySDPSolver

max_forks_repo_name: TomaszWasTaken/MySDPSolver