text:
@testset "fourier_diff" begin
    @test fourier_diff(5, order=1) ≈ [0.0 0.8506508083520398 -0.5257311121191336 0.5257311121191336 -0.8506508083520399; -0.8506508083520399 0.0 0.8506508083520398 -0.5257311121191336 0.5257311121191336; 0.5257311121191336 -0.8506508083520399 0.0 0.8506508083520398 -0.5257311121191336; -0.5257311121191336 0.5257311121191336 -0.8506508083520399 0.0 0.8506508083520398; 0.8506508083520398 -0.5257311121191336 0.5257311121191336 -0.8506508083520399 0.0]
    @test fourier_diff(6, order=1) ≈ [0.0 0.8660254037844387 -0.2886751345948127 -3.061616997868383e-17 0.288675134594813 -0.8660254037844387; -0.8660254037844387 0.0 0.8660254037844387 -0.2886751345948127 -3.061616997868383e-17 0.288675134594813; 0.288675134594813 -0.8660254037844387 0.0 0.8660254037844387 -0.2886751345948127 -3.061616997868383e-17; -3.061616997868383e-17 0.288675134594813 -0.8660254037844387 0.0 0.8660254037844387 -0.2886751345948127; -0.2886751345948127 -3.061616997868383e-17 0.288675134594813 -0.8660254037844387 0.0 0.8660254037844387; 0.8660254037844387 -0.2886751345948127 -3.061616997868383e-17 0.288675134594813 -0.8660254037844387 0.0]
    @test fourier_diff(BigFloat, 5, order=1) ≈ [BigFloat("0.0") BigFloat("0.850650808352039932181540497063011072240401403764816881836740242377884047363955") BigFloat("-0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("-0.8506508083520399321815404970630110722404014037648168818367402423778840473639636"); BigFloat("-0.8506508083520399321815404970630110722404014037648168818367402423778840473639636") BigFloat("0.0") BigFloat("0.850650808352039932181540497063011072240401403764816881836740242377884047363955") BigFloat("-0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("0.5257311121191336060256690848478766072854979322433417815289355232412111464032185"); BigFloat("0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("-0.8506508083520399321815404970630110722404014037648168818367402423778840473639636") BigFloat("0.0") BigFloat("0.850650808352039932181540497063011072240401403764816881836740242377884047363955") BigFloat("-0.5257311121191336060256690848478766072854979322433417815289355232412111464032185"); BigFloat("-0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("-0.8506508083520399321815404970630110722404014037648168818367402423778840473639636") BigFloat("0.0") BigFloat("0.850650808352039932181540497063011072240401403764816881836740242377884047363955"); BigFloat("0.850650808352039932181540497063011072240401403764816881836740242377884047363955") BigFloat("-0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("0.5257311121191336060256690848478766072854979322433417815289355232412111464032185") BigFloat("-0.8506508083520399321815404970630110722404014037648168818367402423778840473639636") BigFloat("0.0")]
    @test fourier_diff(BigFloat, 6, order=1) ≈ [BigFloat("0.0") BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383") BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749") BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78") BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619") BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902"); BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902") BigFloat("0.0") BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383") BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749") BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78") BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619"); BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619") BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902") BigFloat("0.0") BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383") BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749") BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78"); BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78") BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619") BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902") BigFloat("0.0") BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383") BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749"); BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749") BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78") BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619") BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902") BigFloat("0.0") BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383"); BigFloat("0.8660254037844386467637231707529361834714026269051903140279034897259665084543383") BigFloat("-0.2886751345948128822545743902509787278238008756350634380093011632419888361514749") BigFloat("-2.742293602448380191855326565989245052626895591271719877947514292480356721283385e-78") BigFloat("0.2886751345948128822545743902509787278238008756350634380093011632419888361514619") BigFloat("-0.8660254037844386467637231707529361834714026269051903140279034897259665084543902") BigFloat("0.0")]

    @test fourier_diff(5, order=2) ≈ [-2.0 1.1708203932499366 -0.1708203932499369 -0.17082039324993695 1.170820393249937; 1.170820393249937 -2.0 1.1708203932499366 -0.1708203932499369 -0.17082039324993695; -0.17082039324993695 1.170820393249937 -2.0 1.1708203932499366 -0.1708203932499369; -0.1708203932499369 -0.17082039324993695 1.170820393249937 -2.0 1.1708203932499366; 1.1708203932499366 -0.1708203932499369 -0.17082039324993695 1.170820393249937 -2.0]
    @test fourier_diff(6, order=2) ≈ [-3.1666666666666665 2.000000000000001 -0.6666666666666666 0.5 -0.6666666666666669 2.000000000000001; 2.000000000000001 -3.1666666666666665 2.000000000000001 -0.6666666666666666 0.5 -0.6666666666666669; -0.6666666666666669 2.000000000000001 -3.1666666666666665 2.000000000000001 -0.6666666666666666 0.5; 0.5 -0.6666666666666669 2.000000000000001 -3.1666666666666665 2.000000000000001 -0.6666666666666666; -0.6666666666666666 0.5 -0.6666666666666669 2.000000000000001 -3.1666666666666665 2.000000000000001; 2.000000000000001 -0.6666666666666666 0.5 -0.6666666666666669 2.000000000000001 -3.1666666666666665]
    @test fourier_diff(BigFloat, 5, order=2) ≈ [BigFloat("-1.999999999999999999999999999999999999999999999999999999999999999999999999999983") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691329") BigFloat("-0.170820393249936908922752100619382870632185507883457717281269173623156277691348") BigFloat("-0.1708203932499369089227521006193828706321855078834577172812691736231562776913437") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691346"); BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691346") BigFloat("-1.999999999999999999999999999999999999999999999999999999999999999999999999999983") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691329") BigFloat("-0.170820393249936908922752100619382870632185507883457717281269173623156277691348") BigFloat("-0.1708203932499369089227521006193828706321855078834577172812691736231562776913437"); BigFloat("-0.1708203932499369089227521006193828706321855078834577172812691736231562776913437") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691346") BigFloat("-1.999999999999999999999999999999999999999999999999999999999999999999999999999983") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691329") BigFloat("-0.170820393249936908922752100619382870632185507883457717281269173623156277691348"); BigFloat("-0.170820393249936908922752100619382870632185507883457717281269173623156277691348") BigFloat("-0.1708203932499369089227521006193828706321855078834577172812691736231562776913437") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691346") BigFloat("-1.999999999999999999999999999999999999999999999999999999999999999999999999999983") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691329"); BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691329") BigFloat("-0.170820393249936908922752100619382870632185507883457717281269173623156277691348") BigFloat("-0.1708203932499369089227521006193828706321855078834577172812691736231562776913437") BigFloat("1.170820393249936908922752100619382870632185507883457717281269173623156277691346") BigFloat("-1.999999999999999999999999999999999999999999999999999999999999999999999999999983")]
    @test fourier_diff(BigFloat, 6, order=2) ≈ [BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868") BigFloat("0.50") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965"); BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965") BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868") BigFloat("0.50") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609"); BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965") BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868") BigFloat("0.50"); BigFloat("0.50") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965") BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868"); BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868") BigFloat("0.50") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965") BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793"); BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999793") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666868") BigFloat("0.50") BigFloat("-0.6666666666666666666666666666666666666666666666666666666666666666666666666666609") BigFloat("1.999999999999999999999999999999999999999999999999999999999999999999999999999965") BigFloat("-3.166666666666666666666666666666666666666666666666666666666666666666666666666678")]

    @test fourier_diff(5) == fourier_diff(5, order=1)
end


text:
using Luxor

struct SVG <: Drawable
    path::String
    scale::Float16
end

function SVG(
    path;
    scale=1.0,
    )
    SVG(
        path,
        scale,
    )
end

function draw(s::SVG, pos, cell_width; halign=:left, valign=:middle)
    @layer begin
        translate(pos)
        scale(s.scale, s.scale)
        placeimage(
            readsvg(s.path);
            centered=true
        )
    end
end

ext: jl


text:
export hasEqual, hasIdentical, flatten, markUnique, getUnique!

using Statistics: std, mean
using Symbolics
using LinearAlgebra: eigvals, svdvals, eigen

# Function for submodule loading and integrity checking.
function tryIncluding(subModuleName::String; subModulePath=(@__DIR__)[:]*"/SubModule")
    try
        include(subModulePath*"/"*subModuleName*".jl")
        return true
    catch err
        warning = """
        Submodule `$(subModuleName)` failed loading and won't be useable:

            $(err)

        `///magenta///However, this DOES NOT affect the functionality of the main module.`
        """
        printStyledInfo(warning, title="WARNING:\n", titleColor=:light_yellow)
        return false
    end
end

"""

    @compareLength inputArg1 inputArg2 argNames::String... -> length(inputArg1)

A macro that checks whether the lengths of 2 `Arrays`/`Tuples` are equal. It returns the
lengths of the compared objects are the same; it throws a detailed ERROR message when they
are not equal.

You can specify the name for the compared variables in arguments for better ERROR
information.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
julia> Quiqbox.@compareLength [1,2] [3,4]
2

julia> let a = [1,2], b = [3]
           Quiqbox.@compareLength a b
       end
ERROR: The lengths of a and b are NOT equal.
       a::Vector{Int64}   length: 2
       b::Vector{Int64}   length: 1

julia> Quiqbox.@compareLength [1,2] [3] "a"
ERROR: The lengths of a ([1, 2]) and [3] are NOT equal.
       a ([1, 2])::Vector{Int64}   length: 2
       [3]::Vector{Int64}   length: 1

julia> Quiqbox.@compareLength [1,2] [3] "a" "b"
ERROR: The lengths of a ([1, 2]) and b ([3]) are NOT equal.
       a ([1, 2])::Vector{Int64}   length: 2
       b ([3])::Vector{Int64}   length: 1
```
"""
macro compareLength(inputArg1, inputArg2, argNames::String...)
    # In a macro you must escape all the user inputs once and exactly once.
    ns0 = [string(inputArg1), string(inputArg2)]
    ns = ns0
    quote
        local arg1 = $(esc(inputArg1))
        local arg2 = $(esc(inputArg2))
        type = Union{AbstractArray, Tuple}
        (!(arg1 isa type && arg2 isa type)) && error("The compared objects have to be ",
                                                     "Arrays or Tuples!\n")
        if length(arg1) != length(arg2)
            for i = 1:length($argNames)
                # Replace the default type control ERROR message.
                !($argNames[i] isa String) && error("The object's name has to be a "*
                                                    "`String`!\n")
                $ns[i] = $argNames[i]*" ($($ns0[i]))"
            end
            error("""The lengths of $($ns[1]) and $($ns[2]) are NOT equal.
                           $($ns0[1])::$(typeof(arg1))   length: $(length(arg1))
                           $($ns0[2])::$(typeof(arg2))   length: $(length(arg2))
                    """)
        end
        length(arg1)
    end
end

"""

    hasBoolRelation(boolOp::Function, obj1, obj2; ignoreFunction::Bool=false,
                    ignoreContainer::Bool=false, decomposeNumberCollection::Bool=false) ->
    Bool

Recursively apply the specified boolean operator to all the fields within 2 objects
(normally 2 `struct`s in the same type). It returns `true` only if all comparisons
performed return `true`. Note that the boolean operator should have method(s) defined for
all the possible elements inside the compared objects.

If `ignoreFunction = true`, the function will ignore comparisons between Function-type
fields.

If `ignoreContainer = true`, the function will ignore the difference of the container(s) as
long as the boolean operator returns true for the field(s)/entry(s) from two objects
respectively.

If `decomposeNumberCollection = true`, then `Tuple{Vararg{Number}}` and `Array{<:Number}`
will be treated as decomposable containers.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Float64
    end

    a = S(1, 1.0)
    b = S(2, 0.5)
    c = S(2, 1.5)

    Quiqbox.hasBoolRelation(>, a, b) |> println
    Quiqbox.hasBoolRelation(>, b, a) |> println
    Quiqbox.hasBoolRelation(>, c, a) |> println
end

# output
false
false
true
```

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Float64
    end

    struct S2
        a::Int
        b::Float64
    end

    Quiqbox.hasBoolRelation(==, S(1,2), S2(1,2), ignoreContainer=true)
end

# output
true
```
"""
function hasBoolRelation(boolOp::F, obj1, obj2;
                         ignoreFunction::Bool=false,
                         ignoreContainer::Bool=false,
                         decomposeNumberCollection::Bool=false) where {F<:Function}
    res = true
    t1 = typeof(obj1)
    t2 = typeof(obj2)
    if (t1 <: Function) && (t2 <: Function)
        ignoreFunction ? (return true) : (return boolOp(obj1, obj2))
    elseif (t1 <: Number) && (t2 <: Number)
        return boolOp(obj1, obj2)
    elseif t1 != t2 && !ignoreContainer
        return false
    elseif obj1 isa Union{Array, Tuple}
        if !decomposeNumberCollection &&
           (eltype(obj1) <: Number) && (eltype(obj2) <: Number)
            return boolOp(obj1, obj2)
        end
        length(obj1) != length(obj2) && (return false)
        !ignoreContainer && obj1 isa Matrix && (size(obj1) != size(obj2)) && (return false)
        for (i,j) in zip(obj1, obj2)
            res *= hasBoolRelation(boolOp, i, j; ignoreFunction, ignoreContainer,
                                   decomposeNumberCollection)
            !res && (return false)
        end
    elseif obj1 isa Type || obj2 isa Type
        return boolOp(obj1, obj2)
    else
        fs1 = fieldnames(t1)
        fs2 = fieldnames(t2)
        if fs1 == fs2
            if length(fs1) == 0
                res = boolOp(obj1, obj2)
            else
                for i in fs1
                    isdefined(obj1, i) == (fieldDefined = isdefined(obj2, i)) &&
                    (fieldDefined ? nothing : continue)
                    res *= hasBoolRelation(boolOp, getproperty(obj1, i),
                                           getproperty(obj2, i); ignoreFunction,
                                           ignoreContainer, decomposeNumberCollection)
                    !res && (return false)
                end
            end
        else
            return false
        end
    end
    res
end
## Refer overload for `ParamBox` to Overload.jl.

"""

    hasBoolRelation(boolOp::F, obj1, obj2, obj3...;
                    ignoreFunction::Bool=false,
                    ignoreContainer::Bool=false,
                    decomposeNumberCollection::Bool=false) where {F<:Function} ->
    Bool

Method for more than 2 objects. E.g.: `hasBoolRelation(>, a, b, c)` is equivalent to
`hasBoolRelation(>, a, b) && hasBoolRelation(>, b, c)`.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Float64
    end

    a = S(1, 1.0)
    b = S(2, 0.5)
    c = S(2, 1.5)
    d = S(3, 2.0)

    Quiqbox.hasBoolRelation(>=, c, b, a) |> println
    Quiqbox.hasBoolRelation(>=, d, c, b) |> println
end

# output
false
true
```
"""
function hasBoolRelation(boolOp::F, obj1, obj2, obj3...;
                         ignoreFunction::Bool=false,
                         ignoreContainer::Bool=false,
                         decomposeNumberCollection::Bool=false) where {F<:Function}
    res = hasBoolRelation(boolOp, obj1, obj2; ignoreFunction, ignoreContainer)
    tmp = obj2
    if res
        for i in obj3[1:end]
            res *= hasBoolRelation(boolOp, tmp, i; ignoreFunction, ignoreContainer,
                                   decomposeNumberCollection)
            !res && break
            tmp = i
        end
    end
    res
end

"""

    hasEqual(obj1, obj2, obj3...;
             ignoreFunction::Bool=false,
             ignoreContainer::Bool=false,
             decomposeNumberCollection::Bool=false) ->
    Bool

Compare if two objects are the equal.

If `ignoreFunction = true`, the function will ignore comparisons between Function-type
fields.

If `ignoreContainer = true`, the function will ignore the difference of the container(s)
and only compare the field(s)/entry(s) from two objects respectively.

If `decomposeNumberCollection = true`, then `Tuple{Vararg{Number}}` and `Array{<:Number}`
will be treated as decomposable containers.

This is an instantiation of `Quiqbox.hasBoolRelation`.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Float64
    end
    a = S(1, 1.0)
    b = S(1, 1.0)
    c = S(1, 1.0)
    d = S(1, 1.1)

    hasEqual(a, b, c) |> println
    hasEqual(a, b, c, d) |> println
end

# output
true
false
```
"""
hasEqual(obj1, obj2, obj3...; ignoreFunction::Bool=false, ignoreContainer::Bool=false,
        decomposeNumberCollection::Bool=false) =
hasBoolRelation(==, obj1, obj2, obj3...; ignoreFunction, ignoreContainer,
                decomposeNumberCollection)

"""

    hasIdentical(obj1, obj2, obj3...;
                 ignoreFunction::Bool=false,
                 ignoreContainer::Bool=false,
                 decomposeNumberCollection::Bool=false) ->
    Bool

Compare if two objects are the Identical. An instantiation of `hasBoolRelation`.

If `ignoreFunction = true`, the function will ignore comparisons between Function-type
fields.

If `ignoreContainer = true`, the function will ignore the difference of the container(s)
and only compare the field(s)/entry(s) from two objects respectively.

If `decomposeNumberCollection = true`, then `Tuple{Vararg{Number}}` and `Array{<:Number}`
will be treated as decomposable containers.

This is an instantiation of `Quiqbox.hasBoolRelation`.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Array{Float64, 1}
    end

    a = S(1, [1.0, 1.1])
    b = a
    c = b
    d = S(1, [1.0, 1.1])

    hasIdentical(a, b, c) |> println
    hasIdentical(a, b, c, d) |> println
end

# output
true
false
```
"""
hasIdentical(obj1, obj2, obj3...; ignoreFunction::Bool=false, ignoreContainer::Bool=false,
             decomposeNumberCollection::Bool=false) =
hasBoolRelation(===, obj1, obj2, obj3...; ignoreFunction, ignoreContainer,
                decomposeNumberCollection)

"""

    printStyledInfo(str::String; title::String="", titleColor::Symbol=:light_blue) ->
    Nothing

Print info with colorful title and automatically highlighted code blocks enclosed by ` `.

If you want to highlight other contents in different colors, you can also put them inside
` ` and start it with "///theColorSymbolName///". The available color names follows the
values of `color` keyword argument in function `printstyled`.

NOTE: There can only be one color in one ` ` quote.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
julia> Quiqbox.printStyledInfo("This `///magenta///word` is in color magenta.")
This word is in color magenta.
```
"""
function printStyledInfo(str::String;
                         title::String="", titleColor::Symbol=:light_blue)
    codeColor = :cyan
    codeQuotes = findall("`", str)
    l = codeQuotes |> length
    printstyled(title, bold=true, color=titleColor)
    if (l > 0) && (l |> iseven)
        print(str[1:codeQuotes[1][]-1])
        for i = 1:2:l-1
            subStr1 = str[codeQuotes[i][]+1:end]
            blockStr = subStr1[1 : findfirst("`", subStr1)[]-1]
            blockColor = codeColor
            startLoc = findfirst("///", blockStr)
            if startLoc !== nothing
                endLoc = findnext("///", blockStr, 2)
                if endLoc !== nothing
                    colorKey = Symbol(blockStr[startLoc[3]+1 : endLoc[1]-1])
                    blockColor = colorKey
                    blockStr = blockStr[1 : startLoc[1]-1]*blockStr[endLoc[3]+1 : end]
                end
            end
            printstyled(blockStr, color=blockColor)
            subStr2 = str[codeQuotes[i+1][]+1:end]
            loc = findfirst("`", subStr2)
            unBlockStr = subStr2[1 : ((loc===nothing) ? (subStr2|>length) : (loc[]-1))]
            print(unBlockStr)
        end
    else
        print(str)
    end
    println()
end

"""
    flatten(a::Tuple) -> Tuple

    flatten(a::Array) -> Array

Flatten `a::Union{Array, Tuple}` that contains `Array`s and/or `Tuple`s. Only operate on
the outermost layer.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
julia> flatten((:one, 2, [3, 4.0], ([5], "six"), "7"))
(:one, 2, 3.0, 4.0, [5], "six", "7")

julia> flatten([:one, 2, [3, 4.0], ([5], "six"), "7"])
7-element Vector{Any}:
  :one
 2
 3.0
 4.0
  [5]
  "six"
  "7"
```
"""
function flatten(c::Array)
    c2 = map( x->(x isa Union{Array, Tuple} ? x : (x,)), c )
    [(c2...)...]
end

function flatten(c::Tuple)
    c2 = map( x->(x isa Union{Array, Tuple} ? x : (x,)), c )
    ((c2...)...,)
end

"""

    arrayAlloc(arrayLength::Int,
               anExampleOrType::Union{T, Type{T}}) where {T<:Real} ->
    Ptr{T}

Allocate the memory for an array of specified length and element type, then return the
pointer `Ptr` to it.
"""
function arrayAlloc(arrayLength::Int, elementType::Type{T}) where {T<:Real}
    memoryLen = arrayLength*sizeof(elementType) |> Cint
    ccall(:malloc, Ptr{T}, (Cint,), memoryLen)
end

arrayAlloc(arrayLength::Int, NumberExample::T) where {T<:Real} =
arrayAlloc(arrayLength, typeof(NumberExample))

"""

    ArrayPointer{T, N} <: Any

Stores a pointer to the actual address of an array.

≡≡≡ Field(s) ≡≡≡

`ptr::Ptr{T}`: Pointer pointing to the memory address of the first element of the array.

`arr::Array{T, N}`: The mutable array linked to the pointer. As long as the pointer
                    (memory) is not freed, the array is safely preserved.

≡≡≡ Initialization Method(s) ≡≡≡

    ArrayPointer(arr::Array{<:Real, N},
                 showReminder::Bool=true) where {N} -> ArrayPointer{T, N}

Create a `ArrayPointer` that contains a `Ptr` pointing to the actual memory address of the
(1st element of the) `Array`.

To avoid memory leaking, the user should use `free(x.ptr)` after the usage of
`x::ArrayPointer` to free the occupied memory.

If `showReminder=true`, the constructor will pop up a message to remind the user of
such operation.
"""
struct ArrayPointer{T, N} <: Any
    ptr::Ptr{T}
    arr::Array{T, N}

    function ArrayPointer(arr::Array{<:Real, N}, showReminder::Bool=true) where {N}
        len = length(arr)
        elt =  eltype(arr)
        ptr = arrayAlloc(len, elt)
        unsafe_copyto!(ptr, pointer(arr |> copy), len)
        arr2 = unsafe_wrap(Array, ptr, size(arr))
        showReminder && printStyledInfo("""
            Generating a C-array pointer-like object x`::ArrayPointer{$(elt)}`...
            Remember to use free(x.ptr) afterwards to prevent potential memory leaking.
            """)
        new{elt, N}(ptr, arr2)
    end
end

"""

    markUnique(arr::AbstractArray, args...;
               compareFunction::Function = hasEqual, kws...) ->
    markingList:: Array{Int, 1}, uniqueList::Array

Return a `markingList` using `Int` number to mark each different elements from (and inside)
the input argument(s) and a `uniqueList` to contain all the unique elements when
`compareFunction` is set to `hasEqual` (in default).

`args` and `kws` are positional arguments and keywords arguments respectively as
parameters of the specified `compareFunction`.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
markUnique([1, [1, 2],"s", [1, 2]])

# output
([1, 2, 3, 2], Any[1, [1, 2], "s"])
```

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
begin
    struct S
        a::Int
        b::Float64
    end

    a = S(1, 2.0)
    b = S(1, 2.0)
    c = S(1, 2.1)
    d = a

    markUnique([a,b,c,d])
end

# output
([1, 1, 2, 1], S[S(1, 2.0), S(1, 2.1)])
```
"""
function markUnique(arr::AbstractArray, args...;
                    compareFunction::F=hasEqual, kws...) where {F<:Function}
    @assert length(arr) >= 1 "The length of input array should be not less than 1."
    f = (b...)->compareFunction((b..., args...)...; kws...)
    res = Int[1]
    cmprList = eltype(arr)[arr[1]]
    for i = 2:length(arr)
        local j
        isNew = true
        for outer j = 1:length(cmprList)
            if f(cmprList[j], arr[i])
                isNew = false
                break
            end
            j += 1
        end
        push!(res, j)
        isNew && push!(cmprList, arr[i])
    end
    res, cmprList
end

"""

    getUnique!(arr::Array, args...; compareFunction::F = hasEqual, kws...) where
              {F<:Function} ->
    arr::Array

Similar to [`markUnique`](@ref) but instead, just directly return the input `Array` with
repeated entries deleted.

≡≡≡ Example(s) ≡≡≡

```jldoctest; setup = :(push!(LOAD_PATH, "../../src/"); using Quiqbox)
julia> arr = [1, [1, 2],"s", [1, 2]]
4-element Vector{Any}:
 1
  [1, 2]
  "s"
  [1, 2]

julia> getUnique!(arr);

julia> arr
3-element Vector{Any}:
 1
  [1, 2]
  "s"
```
"""
function getUnique!(arr::Array, args...;
                    compareFunction::F = hasEqual, kws...) where {F<:Function}
    @assert length(arr) > 1 "The length of input array should be larger than 1."
    f = (b...)->compareFunction((b..., args...)...; kws...)
    cmprList = eltype(arr)[arr[1]]
    delList = Bool[false]
    for i = 2:length(arr)
        isNew = true
        for j = 1:length(cmprList)
            if f(cmprList[j], arr[i])
                isNew = false
                break
            end
        end
        isNew && push!(cmprList, arr[i])
        push!(delList, !isNew)
    end
    deleteat!(arr, delList)
end

"""
A dummy function that only returns its argument.
"""
itself(x) = x

"""
Similar as `replace` but for Symbols.
"""
function symbolReplace(sym::Symbol, pair::Pair{String, String}; count::Int=typemax(Int))
    replace(sym |> string, pair; count) |> Symbol
end

function renameFunc(fName::String, f::F) where {F<:Function}
    @eval ($(Symbol(fName)))(a...; b...) = $f(a...; b...)
end

"""
Recursively find the final value using the value of each iteration as the key for the
next search.
"""
function recursivelyGet(dict::Dict, startKey::Any)
    res = nothing
    val = get(dict, startKey, missing)
    while !(val isa Missing)
        res = val
        val = get(dict, val, missing)
    end
    res
end

function isOscillateConverged(sequence::Vector{<:Real},
                              threshold1::Real, threshold2::Real=threshold1;
                              leastCycles::Int=1, nPartition::Int=5, returnStd::Bool=false)
    @assert leastCycles>0 && nPartition>1
    len = length(sequence)
    len < leastCycles && (return false)
    slice = len ÷ nPartition
    lastPortion = sequence[max(end-slice, 1) : end]
    remained = sort(lastPortion)[end÷2+1 : end]
    b = std(remained) < threshold1 && abs(sequence[end] - mean(remained)) < threshold2
    returnStd ? (b, std(lastPortion)) : b
end

function splitTerm(term::Symbolics.Num)
    r1 = Symbolics.@rule +(~(~xs)) => [i for i in ~(~xs)]
    r2 = Symbolics.@rule *(~(~xs)) => [[i for i in ~(~xs)] |> prod]
    for r in [r1, r2]
        term = Symbolics.simplify(term, rewriter = r)
    end
    # Converting Symbolics.Arr to Base.Array
    if term isa Symbolics.Arr
        terms = term |> collect
    else
        terms = [term]
    end
    terms
end

function groupedSort(v::Vector, sortFunction::F=itself) where {F<:Function}
    sortedArr = sort(v, by=x->sortFunction(x))
    state1 = 1
    groups = typeof(v)[]
    next = iterate(sortedArr)
    while next !== nothing
        item, state = next
        next = iterate(sortedArr, state)
        if next === nothing || sortFunction(next[1]) != sortFunction(item)
            push!(groups, sortedArr[state1:state-1])
            state1 = state
        end
    end
    groups
end

function mapPermute(arr, permFunction)
    ks = [[true, x, i] for (x, i) in zip(arr, eachindex(arr))]
    arrNew = permFunction(arr)
    idx = Int[]
    for ele in arrNew
        i = findfirst(x -> x[1] == true && hasIdentical(x[2], ele), ks)
        push!(idx, i)
        ks[i][1] = false
    end
    idx
end

function getFunc(fSym::Symbol, failedResult=missing)
    try
        getfield(Quiqbox, fSym)
    catch
        try
            getfield(Main, fSym)
        catch
            try
                fSym |> string |> Meta.parse |> eval
            catch
                (_) -> failedResult
            end
        end
    end
end

struct Pf{C, F} <: ParameterizedFunction
    f::Function
end

Pf(c::Float64, f::Function) = Pf{c, nameOf(f)}(f)
Pf(c::Float64, f::Pf{C, F}) where {C, F} = Pf{c*C, F}(f.f)
Pf(c::Float64, ::Val{T}) where {T} = Pf{c, T}(getFunc(T, NaN))
Pf(c::Float64, ::Val{Pf{C, F}}) where {C, F} = Pf{c*C, F}(getFunc(F, NaN))

(f::Pf{C})(x::Real) where {C} = C * f.f(x)
(::Type{Pf{C, F}})(x::Real) where {C, F} = C * getFunc(F, NaN)(x)

Pf(c::Float64, ::Pf{C, :itself}) where {C} = Pf{c*C, :itself}(itself)
Pf(c::Float64, ::Val{:itself}) = Pf{c, :itself}(itself)
Pf(c::Float64, ::Val{Pf{C, :itself}}) where {C} = Pf{c*C, :itself}(itself)

(f::Pf{C, :itself})(x::Real) where {C} = C * x
(::Type{Pf{C, :itself}})(x::Real) where {C} = C * x

nameOf(f::ParameterizedFunction) = typeof(f)

nameOf(f) = nameof(f)


text:
using Statistics

@testset "Auxiliary Functions Test" begin

    @testset "check constant columns" begin

		try
			PLSRegressor.check_constant_cols([1.0 1;1 2;1 3])
		catch
			@test true
		end
		try
			PLSRegressor.check_constant_cols([1.0;1;1][:,:])
		catch
			@test true
		end
		try
			PLSRegressor.check_constant_cols([1.0 2 3])
		catch
			@test true
		end
		try
			PLSRegressor.check_constant_cols([1.0; 1; 1][:,:])
		catch
			@test true
		end

		@test PLSRegressor.check_constant_cols([1.0 1;2 2;3 3])
		@test PLSRegressor.check_constant_cols([1.0;2;3][:,:])

	end

	@testset "centralize" begin

		X        = [1; 2; 3.0][:,:]
		X        = PLSRegressor.centralize_data(X,mean(X,dims=1),std(X,dims=1))
		@test all(X .== [-1;0;1.0])

	end

	@testset "decentralize" begin

		Xo        = [1; 2; 3.0][:,:]
		Xn        = [-1;0;1.0][:,:]
		Xn        = PLSRegressor.decentralize_data(Xn,mean(Xo,dims=1),std(Xo,dims=1))
		@test all(Xn .== [1; 2; 3.0])

	end

	@testset "checkdata" begin

         try
			 PLSRegressor.check_params(2,1,"linear")
		 catch
			 @test true
		 end
		 try
			 PLSRegressor.check_params(-1,2,"linear")
		 catch
			 @test true
		 end
		 try
			 PLSRegressor.check_params(1,2,"x")
		 catch
			 @test true
		 end

		 @test PLSRegressor.check_params(1,2,"linear")

	end

	@testset "checkparams" begin

		 try
			 PLSRegressor.check_data(zeros(0,0), 0)
		 catch
			 @test true
		 end
		 try
			 PLSRegressor.check_data(zeros(1,1), 10)
		 catch
			 @test true
		 end
		 @test PLSRegressor.check_data(zeros(1,1), 1)

	end

end;


text:
include("../src/FameSVD.jl")

using BenchmarkTools
using LinearAlgebra
using Statistics
using Plots

function benchmark_pd()
  t = Array{Float64, 2}(undef, 10, 2)
  n = 1000

  FameSVD.fsvd(rand(Float64, 4, 2))
  FameSVD.svd(rand(Float64, 4, 2))

  for i = 1:10
    a = @benchmark FameSVD.fsvd(A) setup=(A=randn(Float64, $i * $n, $n))
    b = @benchmark svd(A) setup=(A=randn(Float64, $i * $n, $n))

    t[i, 1] = median(a).time
    t[i, 2] = median(b).time
  end

  png(plot(1:10, t, xaxis="Matrix rows size(x * 1000)", yaxis="execution time", labels=["fsvd", "svd"]), "bench.png")
end


text:
using Yao
using YaoExtensions

# arXiv: 1704.05018, table S2
function hydrogen_hamiltonian()
    Z1 = put(2,1=>Z)
    Z2 = put(2,2=>Z)
    X1 = put(2,1=>X)
    X2 = put(2,2=>X)
    0.011280*Z1*Z2 + 0.397936*Z1 + 0.397936*Z2 + 0.180931*X1*X2
end

using Flux: Optimise
function train!(circ, hamiltonian; optimizer, niter::Int=100)
     params = parameters(circ)
     dispatch!(circ, :random)
     for i=1:niter
         _, grad = expect'(hamiltonian, zero_state(nqubits(circ)) => circ)
         Optimise.update!(optimizer, params, grad)
         dispatch!(circ, params)
         println("Energy = $(expect(hamiltonian, zero_state(nqubits(h)) |> circ) |> real)")
     end
     return expect(hamiltonian, zero_state(nqubits(h)) |> circ)
end

h = hydrogen_hamiltonian()
c = variational_circuit(2)
emin_vqe = train!(c, h; optimizer=Optimise.ADAM(0.1))

using LinearAlgebra
emin = eigvals(Matrix(mat(h)))[1]
@assert isapprox(emin, emin, atol=1e-1)


text:
function calDiagnostics!(
    ocn :: Ocean
)

    ocn.diag[:total_heat][1] = (c_p_sw * ocn.ASUM.sumΔvol_T * reshape(ocn.Ts_mixed, :, 1))[1]
    ocn.diag[:total_salt][1] =          (ocn.ASUM.sumΔvol_T * reshape(ocn.Ss_mixed, :, 1))[1]

end

max_stars_repo_name: meteorologytoday/IOM

max_issues_repo_name: meteorologytoday/IOM

max_forks_repo_name: meteorologytoday/IOM


text:
using OneHotArrays
using Test

@testset "OneHotArray" begin
  include("array.jl")
end

@testset "Constructors" begin
  include("onehot.jl")
end

@testset "Linear Algebra" begin
  include("linalg.jl")
end


text:
using LinearAlgebra

# k_local(Elastic Modulus,thickness,height,length,number of nodes,row/column one,row/column two, ==1 adjusts for fixed node)

#creates local stiffness matrix and maps to global
function k_local(E,t,h,l,i,n1,n2,fixed=0)
    k_global = zeros(i,i)
    A = t*h
    k = ((E*A)/l)
    k_global[n1,n1] = k
    k_global[n1,n2] = -k
    k_global[n2,n1] = -k
    k_global[n2,n2] = k
    if fixed == 1
      adj = zeros(1,i)
      adj[n1] = 1
      k_global[n1,:] = adj
      k_global[:,n1] = adj
    end
    return k_global
end

# constants from problem statement
E = 68.9*10^9
t = 0.006
h13 = 0.08
h2 = 0.02
l = 0.025
l2 = 0.1

# combining all globally mapped matrices
kG = k_local(E,t,h13,l,4,1,2,1)+
     k_local(E,t,h2,l2,4,2,3)+
     k_local(E,t,h13,l,4,3,4)

# Force Vector
f = Float64[0 0 0 1800]
F = transpose(f)

# Displacement vector
u = kG\F
println("Nodal Displacement Vector (meters)")
display(u)

# Strain Function
ϵ(x2,x1,l) = (x2 - x1)/l

#Strain Vector
ϵM = [ϵ(u[2],u[1],l)
      ϵ(u[3],u[2],l2)
      ϵ(u[4],u[3],l)]

# Stress Function
σ(E,ϵ) = E*ϵ

# Stress Vector
σM = [σ(E,ϵM[1])
      σ(E,ϵM[2])
      σ(E,ϵM[3])]
println("Element Stress Vector (Pa)")
display(σM)


text:
### Access to the elements of a stack or sub structures
import Base.isempty

## Testing if a Stack is empty
# CompressedStack
function isempty(stack::CompressedStack)
  bfirst = isempty(stack.first_partial) || isempty(stack.first_partial[1])
  bsecond = isempty(stack.second_partial) || isempty(stack.second_partial[1])
  bcompressed = isnull(stack.compressed)
  bexplicit = isempty(stack.first_explicit) && isempty(stack.second_explicit)
  return bfirst && bsecond && bcompressed && bexplicit
end
function isempty(stack::CompressedStack, lvl::Int)
  if lvl < stack.depth
    return isempty(stack.first_partial[lvl])
  else
    return isempty(stack.first_explicit)
  end
end

# NormalStack
function isempty(stack::NormalStack)
  isempty(stack.data)
end

## Reading the top explicit element of a CompressedStack
function top_explicit(stack::CompressedStack)
  if isempty(stack.first_explicit)
    stack.second_explicit[end].index
  else
    stack.first_explicit[end].index
  end
end

## Reading the top element
# Normal Stacks
function top(stack::NormalStack)
  d = stack.data[end]
  return d.index
end

# Compressed Stacks
function top{T}(block::Block{T})
  return block[end].last
end

function top{T}(stack::CompressedStack{T})
  if !isempty(stack.first_partial) && !isempty(stack.first_partial[1])
    top(stack.first_partial[1])
  elseif !isempty(stack.second_partial) && !isempty(stack.second_partial[1])
    top(stack.second_partial[1])
  elseif !isempty(stack.first_explicit) || !isempty(stack.second_explicit)
    top_explicit(stack)
  else
    return get(stack.compressed).last
  end
end

## Reading bottom element
# Block
function bottom{T}(block::Block{T})
  return block[end].first
end

## Update top elements
# Used in pop.jl
function update_top!{T}(block::Block{T}, index::Int)
  block[end].last = index
end
# Used in push.jl
function update_top!{T}(block::Block{T}, subblock::Int, index::Int)
  block[subblock].last = index
end

function update_signature!{T}(sign::Signature{T}, index::Int)
  sign.last = index
end


text:
function (sim::Simulation)(::typeof(get), sym::Symbol, dir::Direction{D}; reshaped::Bool = true)  where D
	F = getproperty(sim,sym)
	Nw = ncells(sim.grid);
	F_vec = F[((D-1)*Nw +1) : (D*Nw)]
	reshaped ? reshape(F_vec,size(sim.grid)) : F_vec
end

function (sim::Simulation{Dim})(sym::Symbol,dir::Direction{D})  where {Dim,D}

    !(sym === :ϵᵣ && x̂ == dir ) || return sim.ϵᵣ_xx
	!(sym === :ϵᵣ && ŷ == dir ) || return sim.ϵᵣ_yy
	!(sym === :ϵᵣ && ẑ == dir ) || return sim.ϵᵣ_zz

	!(sym === :μᵣ && x̂ == dir ) || return sim.μᵣ_xx
	!(sym === :μᵣ && ŷ == dir ) || return sim.μᵣ_yy
	!(sym === :μᵣ && ẑ == dir ) || return sim.μᵣ_zz

	!(sym === :J && x̂ == dir ) || return sim.J_x
	!(sym === :J && ŷ == dir ) || return sim.J_y
	!(sym === :J && ẑ == dir ) || return sim.J_z

	!(sym === :PML && x̂ == dir ) || return sim.S_x
	!(sym === :PML && ŷ == dir ) || return sim.S_y
	!(sym === :PML && ẑ == dir ) || return sim.S_z

	!(sym === :Q && x̂ == dir ) || return sim.Q
	!(sym === :Q && ŷ == dir ) || return sim.Q
	!(sym === :Q && ẑ == dir ) || return sim.Q

	# Electric Field
	!(sym === :E && x̂ == dir ) || return sim(get,:E,x̂)
	!(sym === :E && ŷ == dir ) || return sim(get,:E,ŷ)
	!(sym === :E && ẑ == dir ) || return sim(get,:E,ẑ)

	# Magnetic Field
	!(sym === :H && x̂ == dir ) || return  -1 / (im*2pi/sim.λ₀ ) .* (∂(p̂,sim.grid,ŷ,sim(:E,ẑ)) - ∂(p̂,sim.grid,ẑ,sim(:E,ŷ)))
	!(sym === :H && ŷ == dir ) || return   1 / (im*2pi/sim.λ₀ ) .* (∂(p̂,sim.grid,x̂,sim(:E,ẑ)) - ∂(p̂,sim.grid,ẑ,sim(:E,x̂)))
	!(sym === :H && ẑ == dir ) || return  -1 / (im*2pi/sim.λ₀ ) .* (∂(p̂,sim.grid,x̂,sim(:E,ŷ)) - ∂(p̂,sim.grid,ŷ,sim(:E,x̂)))

	# Complex Poynting Vector
	!(sym === :S && x̂ == dir ) || return  sim(:E,ŷ)  .* conj.(sim(:H,ẑ)) - sim(:E,ẑ) .* conj.(sim(:H,ŷ))
	!(sym === :S && ŷ == dir ) || return -(sim(:E,x̂) .* conj.(sim(:H,ẑ)) - sim(:E,ẑ) .* conj.(sim(:H,x̂)))
	!(sym === :S && ẑ == dir ) || return  sim(:E,x̂)  .* conj.(sim(:H,ŷ)) - sim(:E,ŷ) .* conj.(sim(:H,x̂))
end

function (sim::Simulation{Dim})(sym::Symbol,dir::Direction,r::Tuple,gt::GridType = p̂)  where {Dim}
	inds = CartesianIndices(gt,sim.grid,r)
	return  sim(sym,dir)[inds]
end


text:
if VERSION >= v"0.7"
  if !(eval(:(@isdefined third_float)))
    const third_float = 1.0/3.0
    const third_big= big(1)/big(3)
    const twothird_float = 2.0/3.0
    const twothird_big = big(2)/big(3)
  end
else
  if !(isdefined(:third_float))
    const third_float = 1.0/3.0
    const third_big= big(1)/big(3)
    const twothird_float = 2.0/3.0
    const twothird_big = big(2)/big(3)
  end
end

max_stars_repo_name: tagordon/limbdark

max_issues_repo_name: tagordon/limbdark

max_forks_repo_name: tagordon/limbdark


text:
"""
    sdataset(datasetname)
Load popular spatial dataset. See the package documentation for a list of available datasets.

# Examples
```jldoctest
julia> guerry = sdataset("Guerry");
```
"""
function sdataset(datasetname::AbstractString)::DataFrame

    sdatasetable = nothing
    if datasetname == "Guerry"
        sdatasetable = Shapefile.Table(joinpath(artifact"GeoDaGuerry", "guerry.shp"))
    elseif datasetname == "Bostonhsg"
        sdatasetable = Shapefile.Table(joinpath(artifact"GeoDaBostonhsg", "boston.shp"))
    else
        throw(ArgumentError("Dataset $(datasetname) not available. See the package documentation for a list of available datasets."))
    end

    DataFrame(sdatasetable)
end


text:
# Licensed to the Apache Software Foundation (ASF) under one
# or more contributor license agreements.  See the NOTICE file
# distributed with this work for additional information
# regarding copyright ownership.  The ASF licenses this file
# to you 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.

struct CompressedBuffer
    data::Vector{UInt8}
    uncompressedlength::Int64
end

"""
    Arrow.Compressed

Represents the compressed version of an [`ArrowVector`](@ref).
Holds a reference to the original column. May have `Compressed`
children for nested array types.
"""
struct Compressed{Z, A}
    data::A
    buffers::Vector{CompressedBuffer}
    len::Int64
    nullcount::Int64
    children::Vector{Compressed}
end

Base.length(c::Compressed) = c.len
Base.eltype(c::Compressed{Z, A}) where {Z, A} = eltype(A)
getmetadata(x::Compressed) = getmetadata(x.data)
compressiontype(c::Compressed{Z}) where {Z} = Z

function compress(Z::Meta.CompressionType, comp, x::Array)
    GC.@preserve x begin
        y = unsafe_wrap(Array, convert(Ptr{UInt8}, pointer(x)), sizeof(x))
        return CompressedBuffer(transcode(comp, y), length(y))
    end
end

compress(Z::Meta.CompressionType, comp, x) = compress(Z, comp, convert(Array, x))

compress(Z::Meta.CompressionType, comp, v::ValidityBitmap) =
    v.nc == 0 ? CompressedBuffer(UInt8[], 0) : compress(Z, comp, view(v.bytes, v.pos:(v.pos + cld(v.ℓ, 8) - 1)))

function makenodesbuffers!(col::Compressed, fieldnodes, fieldbuffers, bufferoffset, alignment)
    push!(fieldnodes, FieldNode(col.len, col.nullcount))
    @debug 1 "made field node: nodeidx = $(length(fieldnodes)), col = $(typeof(col)), len = $(fieldnodes[end].length), nc = $(fieldnodes[end].null_count)"
    for buffer in col.buffers
        blen = length(buffer.data) == 0 ? 0 : 8 + length(buffer.data)
        push!(fieldbuffers, Buffer(bufferoffset, blen))
        @debug 1 "made field buffer: bufferidx = $(length(fieldbuffers)), offset = $(fieldbuffers[end].offset), len = $(fieldbuffers[end].length), padded = $(padding(fieldbuffers[end].length, alignment))"
        bufferoffset += padding(blen, alignment)
    end
    for child in col.children
        bufferoffset = makenodesbuffers!(child, fieldnodes, fieldbuffers, bufferoffset, alignment)
    end
    return bufferoffset
end

function writearray(io, b::CompressedBuffer)
    if length(b.data) > 0
        n = Base.write(io, b.uncompressedlength)
        @debug 1 "writing compressed buffer: uncompressedlength = $(b.uncompressedlength), n = $(length(b.data))"
        @debug 2 b.data
        return n + Base.write(io, b.data)
    end
    return 0
end

function writebuffer(io, col::Compressed, alignment)
    @debug 1 "writebuffer: col = $(typeof(col))"
    @debug 2 col
    for buffer in col.buffers
        n = writearray(io, buffer)
        writezeros(io, paddinglength(n, alignment))
    end
    for child in col.children
        writebuffer(io, child, alignment)
    end
    return
end

max_forks_repo_name: XpressAI/arrow


text:
function used_memory()
    int(Sys.total_memory()) - int(Sys.free_memory())
end

function benchmark(count, warmup, target_time, func, gc_before)
    force_gc()
    first_execution = @elapsed func()
end

function run_benchmark(env::Environment,
                       bench::Benchmark,
                       force_gc::Bool)
    # run for at least a 1/10 of a second
    run_for_atleast(0.1, 10000, time_clock)
    min_time = min(env.clock_resolution * 1000, 0.1)

    test_time, test_iter, _ = run_for_atleast(min_time, 1, bench.run)
    @printf("ran %d iterations in %s\n", test_iter, time_str(test_time))

    new_iters = int(ceil(min_time * test_iter / test_time))

    #sample_count <- config
    sample_count = 1000

    gc_time = @elapsed Base.gc()

    est_time = sample_count * new_iters * test_time * bench.ntimes / test_iter
    est_time = force_gc ? est_time + (new_iters * sample_count  * gc_time) : est_time
    if true || est_time > 5
    	@printf("collecting %d samples, %d iterations each, estimated time %s\n",
	    	sample_count, new_iters, time_str(est_time))
    end

    function run_once()
        bench.run(1)
    end

    @printf("\rProgress 0.0%%")
    run_once()
    times = zeros(sample_count)
    time_begin = time_sec()
    for sample in 1:sample_count
        progress = int((sample / sample_count) * 100)
        if progress % 10 == 0
 		time_now = time_sec()
    		diff = time_now - time_begin
		to_go = (diff * (100.0 / progress)) - diff
		if to_go > 1
                    @printf("\rprogress %d%%  |  time remaining %s",
                            progress, time_str(to_go))
                end
        end
        t = timed_noresult(run_once, int(new_iters))
        times[sample] = (t - env.clock_cost) / new_iters
        if force_gc
            Base.gc()
        end
    end
    @printf("\rprogress %d%%  |  time remaining %s\n",
	    100, "0.00000 s")
    return times
end

function report_benchmark(desc::String, est::Estimate)
    @printf("%s: %s, lb %s, ub %s, ci %.3f\n",
            desc, time_str(est.point), time_str(est.lbound),
            time_str(est.ubound), est.confidence_level)
end

function report_outlier_variance(outvar::OutlierVariance)
    effect_str = (e::OutlierEffect) ->
	if     typeof(e) == Severe     "severely inflated"
	elseif typeof(e) == Moderate   "moderately inflated"
        elseif typeof(e) == Slight     "slightly inflated"
        elseif typeof(e) == Unaffected "unaffected"
        end
    @printf("variance introduced by outliers: %.3f%%\n", outvar.frac)
    @printf("variance is %s by outliers\n", effect_str(outvar.effect))
end

function run_analyze_one(env::Environment, desc::String, bench::Benchmark)
    times = run_benchmark(env, bench, true)
    # ci <- get config item
    # num_resamples <- get config item
    tail_quantile = 0.025
    n_resamples = 1000
    @printf("analyzing with %d resamples...\n", n_resamples)
    analysis = analyze_sample(times, tail_quantile, n_resamples)
    report_benchmark("mean ", analysis.mean)
    report_benchmark("stdev", analysis.std)
    outliers = classify_outliers(times)
    note_outliers(outliers)
    report_outlier_variance(analysis.outlier_variance)
    (times, analysis, outliers)
end


text:
using RCall, MixedModels, Test
using StatsBase: zscore
using StatsModels: SeqDiffCoding
using Tables: columntable

import JellyMe4: _set_lmer, _set_afex_installed

const LMM = LinearMixedModel
const GLMM = GeneralizedLinearMixedModel

@testset "lmerMod" begin
    reval("""
    if(!("lme4" %in% installed.packages())){
        # use tmp for tests if lme4 isn't available
        .libPaths("/tmp")
        lib <- .libPaths()[1L]
        install.packages("lme4",repos="https://cloud.r-project.org", libs=lib)
    }
    if(!("afex" %in% installed.packages())){
        .libPaths("/tmp")
        lib <- .libPaths()[1L]
        install.packages("afex",repos="https://cloud.r-project.org", libs=lib)
    }
    """)

    # this is available in MixedModels.dataset(:sleepstudy) but with different
    # capitalization than in R
    sleepstudy = rcopy(R"sleepstudy")

    kb07 = dataset(:kb07)

    @testset "get lmerMod" begin
        ### from R ###

        jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1 + Days|Subject)),sleepstudy), REML=false)
        rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1 + Days|Subject),sleepstudy,REML=FALSE)")

        @test jlmm.θ ≈ rlmm.θ atol=0.001
        @test objective(jlmm) ≈ objective(rlmm) atol=0.001
        @test fixef(jlmm) ≈ fixef(rlmm) atol=0.001

        jlmm = refit!(jlmm, REML=true)
        rlmm = rcopy(R"update(m, REML=TRUE)")

        @test jlmm.θ ≈ rlmm.θ atol=0.001
        @test objective(jlmm) ≈ objective(rlmm) atol=0.001
        @test fixef(jlmm) ≈ fixef(rlmm) atol=0.001

        @testset "merModLmerTest" begin
            jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1 + Days|Subject)),sleepstudy), REML=false)
            rlmm = rcopy(R"m <- lmerTest::lmer(Reaction ~ 1 + Days + (1 + Days|Subject),sleepstudy,REML=FALSE)")

            @test jlmm.θ ≈ rlmm.θ atol=0.001
            @test objective(jlmm) ≈ objective(rlmm) atol=0.001
            @test fixef(jlmm) ≈ fixef(rlmm) atol=0.001

            jlmm = refit!(jlmm, REML=true)
            rlmm = rcopy(R"update(m, REML=TRUE)")

            @test jlmm.θ ≈ rlmm.θ atol=0.001
            @test objective(jlmm) ≈ objective(rlmm) atol=0.001
            @test fixef(jlmm) ≈ fixef(rlmm) atol=0.001
        end

        # single scalar RE is different because RCall converts the single-entry θ
        # to a scalar value, which we have to correct
        # of course, this isn't a problem in the other direction, because
        # the scalar-vector distinction for θ is missing in R
        @testset "scalar RE" begin
            jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + (1|Subject)),sleepstudy), REML=false)
            rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1|Subject),sleepstudy,REML=FALSE)")

            @test jlmm.θ ≈ rlmm.θ atol=0.001
            @test objective(jlmm) ≈ objective(rlmm) atol=0.001
            @test fixef(jlmm) ≈ fixef(rlmm) atol=0.001
        end

        @testset "sorting by n BLUPs vs. n groups" begin
            @rput kb07
            rlmm = rcopy(R"rlmm <- lmer(rt_trunc ~ 1 + (1|subj)+(1+load+prec+spkr|item), $(kb07), REML=FALSE)")
            @test rcopy(R"fitted(rlmm)") ≈ fitted(rlmm)
            @test rcopy(R"deviance(rlmm)") ≈ objective(rlmm)
            @test all(isapprox.(collect(VarCorr(rlmm).σρ.item.σ),
                                rcopy(R"""attr(VarCorr(rlmm)[["item"]], "stddev")""")))
            @test all(isapprox.(collect(VarCorr(rlmm).σρ.subj.σ),
                                rcopy(R"""attr(VarCorr(rlmm)[["subj"]], "stddev")""")))
        end

        @testset "contrasts" begin
            reval("""
            data(cake)
            cake\$rr <- with(cake, replicate:recipe)
            """);
            rlmm = rcopy(R"fm1 <- lme4::lmer(angle ~ recipe * temperature + (1|rr), cake, REML= FALSE)");
            @test fixef(rlmm) ≈ rcopy(R"fixef(fm1)");
            # rlmm = rcopy(R"""fm1 <- lme4::lmer(angle ~ recipe * temperature + (1|rr), cake, REML= FALSE,
            #                              contrasts=list(temperature=contr.helmert))""");
            # @test fixef(rlmm) ≈  rcopy(R"fixef(fm1)");
        end

        @testset "caret" begin
            reval("""
            data(cake)
            cake\$rr <- with(cake, replicate:recipe)
            """);
            rlmm = rcopy(R"fm1 <- lme4::lmer(angle ~ (recipe + temperature)^2 + (1|rr), cake, REML= FALSE)");
            @test fixef(rlmm) ≈ rcopy(R"fixef(fm1)");
        end

        @testset "double-bar" begin
            # double bar works because lme4 treats this internally as a rewrite rule into distinct terms
            # the distinct terms are handled the same way in both languages, so everything is fine.
            # printing differs a lot between languages, but that's life

            reval("""
            machines <- as.data.frame(nlme::Machines)
            mach <- lme4::lmer(score ~ Machine + (Machine || Worker), machines, REML=FALSE)
            """)
            machines = rcopy(R"machines")
            @test_throws ArgumentError rlmm = rcopy(R"mach")

            # note that is the lme4 definition of double || -- it really is just a convenience wrapper for
            # splitting terms up this way!
            # jlmm = fit(MixedModel, @formula(score ~ 1 +  Machine + (1|Worker) + (fulldummy(Machine)|Worker)), machines)
            # as a cheat for comparing the covariance matrices, we use PCA
            # @test only(rlmm.rePCA) ≈ only(jlmm.rePCA) atol=0.05

            rlmm = rcopy(R"m <- lme4::lmer(Reaction ~ 1 + Days + (1 + Days||Subject),sleepstudy,REML=FALSE)")
            jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + zerocorr(1 + Days|Subject)),sleepstudy), REML=false)
            # as a cheat for comparing the covariance matrices, we use PCA
            @test only(rlmm.rePCA) ≈ only(jlmm.rePCA) atol=0.05
        end

        @testset "dummy" begin
            # TODO
        end
    end

    @testset "put lmerMod" begin
        ### from Julia ###
        jlmm = LMM(@formula(Reaction ~ 1 + Days + (1 + Days|Subject)),sleepstudy)
        jm = (jlmm, sleepstudy)
        # unfitted model
        @test_throws ArgumentError @rput jm
        fit!(jlmm, REML=true)
        @rput jm
        @test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
        @test rcopy(R"REMLcrit(jm)") ≈ objective(jlmm)

        refit!(jlmm, REML=false)
        jm = (jlmm, sleepstudy)
        @rput jm
        @test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
        @test rcopy(R"deviance(jm)") ≈ objective(jlmm)

        @testset "columntable" begin
            jm = (jlmm, columntable(sleepstudy));
            @rput jm;
        end

        @testset "transformations" begin
            sleepstudy[!,:Days2] = sleepstudy.Days .+ 1
            @rput sleepstudy;
            R"m <- lme4::lmer(log10(Reaction) ~ 1 + log(Days2) + (1 + log(Days2)|Subject),sleepstudy,REML=FALSE)"
            jlmm = fit!(LMM(@formula(log10(Reaction) ~ 1 + log(Days2) + (1 + log(Days2)|Subject)),sleepstudy), REML=false)
            jm = (jlmm, sleepstudy)
            @rput jm;
            @test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
            @test rcopy(R"deviance(jm)") ≈ objective(jlmm)

            jlmm = fit!(LMM(@formula(Reaction ~ 1 + round(Days) + (1|Subject)),sleepstudy), REML=false)
            jm = (jlmm, sleepstudy);

            @test_throws ArgumentError (@rput jm)
        end

        @testset "sorting by n BLUPs vs. n groups" begin
            jlmm = fit(LMM, @formula(rt_trunc ~ 1 + (1|subj)+(1+load+prec+spkr|item)), kb07)
            jm = (jlmm, kb07)
            @rput jm
            @test rcopy(R"fitted(jm)") ≈ fitted(jlmm)
            @test rcopy(R"deviance(jm)") ≈ objective(jlmm)
            @test all(isapprox.(collect(VarCorr(jlmm).σρ.item.σ),
                                rcopy(R"""attr(VarCorr(jm)[["item"]], "stddev")""")))
            @test all(isapprox.(collect(VarCorr(jlmm).σρ.subj.σ),
                                rcopy(R"""attr(VarCorr(jm)[["subj"]], "stddev")""")))
        end

        @testset "contrasts" begin
            cake = rcopy(reval("""
            data(cake)
            cake\$rr <- with(cake, replicate:recipe)
            cake
            """));
            jlmm = fit(MixedModel, @formula(angle ~ recipe * temperature + (1|rr)),
                       cake, REML=false, contrasts=Dict(:temperature => SeqDiffCoding()));
            jm = (jlmm, cake);
            @rput jm;
            @test fixef(jlmm) ≈  rcopy(R"fixef(jm)");
        end

        @testset "fulldummy" begin
            machines = rcopy(R"as.data.frame(nlme::Machines)")
            jlmm = fit(MixedModel, @formula(score ~ 1 +  Machine + (1 + fulldummy(Machine)|Worker)), machines)
            rlmm = (jlmm, machines)
            @rput rlmm
            rlmmrepca = rcopy(R"summary(rePCA(rlmm))$Worker$importance[3,]")
            @test  rlmmrepca ≈ only(MixedModels.rePCA(jlmm; corr=false)) atol=0.05
        end

        @testset "zerocorr" begin
            # TODO: test fulldummy within a zerocorr when zerocorr is better supported

            @testset "lme4" begin
                _set_lmer("lme4::lmer")
                _set_afex_installed(false)

                jlmm = fit!(LMM(@formula(Reaction ~ 1 + Days + zerocorr(1 + Days|Subject)),sleepstudy), REML=false)
                rlmm = (jlmm, sleepstudy)
                @rput rlmm
                @test rcopy(R"""!is(rlmm,"merModLmerTest")""")
                @test only(ranef(jlmm))' ≈ Matrix(rcopy(R"ranef(rlmm)$Subject"))
                @test fixef(jlmm) ≈ rcopy(R"fixef(rlmm)")
                @test vcov(jlmm) ≈ rcopy(R"as.matrix(vcov(rlmm))")
            end

            @testset "afex" begin

                machines = rcopy(R"as.data.frame(nlme::Machines)")

                jlmm = fit(MixedModel, @formula(score ~ 1 +  Machine + zerocorr(0+Machine|Worker)), machines)

                @testset "afex pre-enabled" begin
                    _set_lmer("afex::lmer_alt")
                    rlmm = (jlmm, machines)
                    @rput rlmm
                    @test only(ranef(jlmm))' ≈ Matrix(rcopy(R"ranef(rlmm)$Worker"))
                    @test fixef(jlmm) ≈ rcopy(R"fixef(rlmm)")
                    @test vcov(jlmm) ≈ rcopy(R"as.matrix(vcov(rlmm))")
                end

                @testset "afex auto-enabled" begin
                    _set_lmer("lme4::lmer")
                    _set_afex_installed(true)
                    rlmm = (jlmm, machines)
                    # match_mode needs to specified in case there are further
                    # warnings from R/RCall
                    @test_logs (:info, r"afex::lmer_alt") match_mode=:any @rput rlmm
                    @test only(ranef(jlmm))' ≈ Matrix(rcopy(R"ranef(rlmm)$Worker"))
                    @test fixef(jlmm) ≈ rcopy(R"fixef(rlmm)")
                    @test vcov(jlmm) ≈ rcopy(R"as.matrix(vcov(rlmm))")
                end

                @testset "afex unavailable" begin
                    _set_lmer("lme4::lmer")
                    _set_afex_installed(false)
                    rlmm = (jlmm, machines)
                    @test_throws ArgumentError @rput rlmm
                end
            end

        end
    end
end


text:
abstract type RackettTranslationModel <: TranslationModel end

struct RackettTranslationParam <: EoSParam
    Vc::SingleParam{Float64}
end

@newmodelsimple RackettTranslation RackettTranslationModel RackettTranslationParam

export RackettTranslation
function RackettTranslation(components::Vector{String}; userlocations::Vector{String}=String[], verbose::Bool=false)
    params = getparams(components, ["properties/critical.csv"]; userlocations=userlocations, verbose=verbose)
    Vc = params["vc"]
    packagedparams = RackettTranslationParam(Vc)
    model = RackettTranslation(packagedparams, verbose=verbose)
    return model
end

function translation(model::CubicModel,V,T,z,translation_model::RackettTranslation)
    Tc = model.params.Tc.values
    Pc = model.params.Pc.values
    Vc = translation_model.params.Vc.values
    c = zeros(eltype(Tc),length(Tc))
    for i ∈ @comps
        Tci = Tc[i]
        Pci = Pc[i]
        RT = Tci*R̄
        Zc = Pci*Vc[i]/RT
        c[i] = 0.40768*RT/Pci*(0.29441-Zc)
    end
    return c
end


text:
using Gadfly

set_default_plot_size(6inch, 3inch)

x, y, size = rand(10), rand(10), rand(["foo","bar","pooh"], 10)
p1 = plot(x=x, y=y, size=size, Geom.point)
p2 = plot(x=x, y=y, size=size, Geom.point, Scale.size_discrete2)

hstack(p1, p2)


text:
@add_benchmark function affine_negate_atom(handle_problem)
    x = Variable()
    p = minimize(-x, [x <= 0])
    handle_problem(p)
end

@add_benchmark function affine_kron_atom(handle_problem)
    x = ComplexVariable(3, 3)
    y = [1.0 2.0; 3.0 4.0]
    p = satisfy(kron(x,y) == kron(eye(3), y))
    handle_problem(p)
end

@add_benchmark function affine_mult_atom(handle_problem)
    probs=[
        let
            x = Variable(1)
            p = minimize(2.0 * x, [x >= 2, x <= 4])
        end,
        let
            x = Variable(2)
            A = 1.5 * eye(2)
            p = minimize([2 2] * x, [A * x >= [1.1; 1.1]])
        end,
        let
            y = Variable(1)
            x = Variable(3)
            z = [1.0, 2.0, 3.0] * y
            k = -y * [1.0, 2.0, 3.0]
            c = [y <= 3.0, y >= 0.0, x >= ones(3),  k <= x, x <= z]
            o = 3 * y
            p = Problem(:minimize, o, c)
        end,
        let
            x = ComplexVariable(2,2)
            p = minimize(   real( [1.0im, 0.0]' * x * [1.0im, 0.0] ),
                    [ x == [1.0 0.0; 0.0 1.0] ] )
        end
    ]
    handle_problem.(probs)
end

@add_benchmark function affine_dot_atom(handle_problem)
    probs = [
        let
            x = Variable(2)
            p = minimize(dot([2.0; 2.0], x), x >= [1.1; 1.1])
        end,
        let
            x = Variable(2,2)
            p = minimize(dot(fill(2.0, (2,2)), x), x >= 1.1)
        end
    ]

    handle_problem.(probs)
end

@add_benchmark function affine_add_atom(handle_problem)

    probs = [
        let
            x = Variable(1)
            y = Variable(1)
            p = minimize(x + y, [x >= 3, y >= 2])
        end,
        let
            x = Variable(1)
            p = minimize(x, [eye(2) + x >= eye(2)])
        end,
        let
            y = Variable()
            p = minimize(y - 5, y >= -1)
        end ]

    handle_problem.(probs)
end

@add_benchmark function affine_transpose_atom(handle_problem)
    probs = [
        let
            x = Variable(2)
            c = ones(2, 1)
            p = minimize(x' * c, x >= 1)
        end,
        let
            X = Variable(2, 2)
            c = ones(2, 1)
            p = minimize(c' * X' * c, [X >= ones(2, 2)])
        end,
        let
            rows = 2
            cols = 3
            r = rand(rows, cols)
            r_2 = rand(cols, rows)
            x = Variable(rows, cols)
            c = ones(1, cols)
            d = ones(rows, 1)
            p = minimize(c * x' * d + d' * x * c' + (c * x''''' * d)',
                        [x' >= r_2, x >= r, x''' >= r_2, x'' >= r])
        end ]

        handle_problem.(probs)
end

@add_benchmark function affine_index_atom(handle_problem)
        probs = [
        let
            x = Variable(2)
            p = minimize(x[1] + x[2], [x >= 1])
        end,
        let
            x = Variable(3)
            I = [true true false]
            p = minimize(sum(x[I]), [x >= 1])
        end,
        let
            rows = 6
            cols = 8
            n = 2
            X = Variable(rows, cols)
            A = randn(rows, cols)
            c = rand(1, n)
            p = minimize(c * X[1:n, 5:5+n-1]' * c', X >= A)
        end]
        handle_problem.(probs)
end

@add_benchmark function affine_sum_atom(handle_problem)
        probs = [
        let
            x = Variable(2,2)
            p = minimize(sum(x) - 2*x[1,1], x>=1, x[1,1]<=2)
        end,
        let
            x = Variable(2,2)
            p = minimize(sum(x) - 2*x[1,1], x>=1, x[1,1]<=2)
        end
        ]
    return handle_problem.(probs)
end

@add_benchmark function affine_diag_atom(handle_problem)
    probs = [
        let
            x = Variable(2,2)
            p = minimize(sum(diag(x,1)), x >= 1)
        end,
        let
            x = Variable(4, 4)
            p = minimize(sum(diag(x)), x >= 2)
        end,
    ]
    return handle_problem.(probs)

end

@add_benchmark function affine_tr_atom(handle_problem)
    x = Variable(2,2)
    p = minimize(tr(x), x >= 1)
    handle_problem(p)
end


text:
function EventBasedManeuverTriggers(arg0::AbstractDetector, arg1::AbstractDetector)
    return EventBasedManeuverTriggers((AbstractDetector, AbstractDetector), arg0, arg1)
end

function event_occurred(obj::EventBasedManeuverTriggers, arg0::SpacecraftState, arg1::EventDetector, arg2::jboolean)
    return jcall(obj, "eventOccurred", Action, (SpacecraftState, EventDetector, jboolean), arg0, arg1, arg2)
end

function get_events_detectors(obj::EventBasedManeuverTriggers)
    return jcall(obj, "getEventsDetectors", Stream, ())
end

function get_field_events_detectors(obj::EventBasedManeuverTriggers, arg0::Field)
    return jcall(obj, "getFieldEventsDetectors", Stream, (Field,), arg0)
end

function get_start_firing_detector(obj::EventBasedManeuverTriggers)
    return jcall(obj, "getStartFiringDetector", AbstractDetector, ())
end

function get_stop_firing_detector(obj::EventBasedManeuverTriggers)
    return jcall(obj, "getStopFiringDetector", AbstractDetector, ())
end

function get_triggered_end(obj::EventBasedManeuverTriggers)
    return jcall(obj, "getTriggeredEnd", AbsoluteDate, ())
end

function get_triggered_start(obj::EventBasedManeuverTriggers)
    return jcall(obj, "getTriggeredStart", AbsoluteDate, ())
end

function init(obj::EventBasedManeuverTriggers, arg0::SpacecraftState, arg1::AbsoluteDate)
    return jcall(obj, "init", void, (SpacecraftState, AbsoluteDate), arg0, arg1)
end

function is_firing(obj::EventBasedManeuverTriggers, arg0::AbsoluteDate)
    return jcall(obj, "isFiring", jboolean, (AbsoluteDate,), arg0)
end

function is_firing(obj::EventBasedManeuverTriggers, arg0::AbsoluteDate, arg1::Vector{jdouble})
    return jcall(obj, "isFiring", jboolean, (AbsoluteDate, Vector{jdouble}), arg0, arg1)
end

function is_firing(obj::EventBasedManeuverTriggers, arg0::FieldAbsoluteDate, arg1::Vector{RealFieldElement})
    return jcall(obj, "isFiring", jboolean, (FieldAbsoluteDate, Vector{RealFieldElement}), arg0, arg1)
end

function reset_state(obj::EventHandler, arg0::EventDetector, arg1::SpacecraftState)
    return jcall(obj, "resetState", SpacecraftState, (EventDetector, SpacecraftState), arg0, arg1)
end

function set_firing(obj::EventBasedManeuverTriggers, arg0::jboolean, arg1::AbsoluteDate)
    return jcall(obj, "setFiring", void, (jboolean, AbsoluteDate), arg0, arg1)
end


text:
using ITensors

#
# DMRG calculation of the extended Hubbard model
# ground state wavefunction, and spin densities
#

# Usage:
#
# First install `ArgParse.jl`:
#
# julia>] add ArgParse
#
# Then run:
#
# $ julia exthubbard.jl input.jl
#
# or with options, like
#
# $ julia exthubbard.jl input.jl --N=100 --Npart=20 --nsweep=4 --cutoff=1e-11 --maxdim=20 40 80 100 --noise=1e-12 1e-13
#
# In the REPL (if you just type `julia`), you can call it as follows:
#
# julia> push!(ARGS, "input.jl");
#
# julia> push!(ARGS, "--N=100", "--Npart=20");
#
# julia> include("exthubbard.jl");
#

# Include the specified input file
# from the first input.
# Searches for the default file input.jl
# if nothing is specified.
filename = get(ARGS, 1, "input.jl")
if isfile(filename)
  include(filename)
else
  include("input.jl")
end

#
# Parse the arguments and store
# them in a dictionary.
# Settings are defined in the input file.
#

args = parse_args(settings)

#
# Get the values from the dictionary
#

N = args["N"]
Npart = args["Npart"]
t1 = args["t1"]
t2 = args["t2"]
U = args["U"]
V1 = args["V1"]
nsweep = args["nsweep"]
maxdim = args["maxdim"]
mindim = args["mindim"]
cutoff = args["cutoff"]
noise = args["noise"]

#
# Alternatively, define all of the
# variables for all of the inputs
# programatically. You can use the
# keyword argument `as_symbols = true`
# in the `parse_args` function to
# automatically turn the keys of
# the dictionary into symbols.
#

#for (arg, val) in args
#  arg == "input_file" && continue
#  @eval $(Symbol(arg)) = $val
#end

sweeps = Sweeps(nsweep)
maxdim!(sweeps, maxdim...)
mindim!(sweeps, mindim...)
cutoff!(sweeps, cutoff...)
noise!(sweeps, noise...)

@show sweeps
@show N, Npart

sites = siteinds("Electron", N; conserve_qns=true)

ampo = AutoMPO()
for i=1:N
  ampo .+= U,"Nupdn",i
end
for b=1:N-1
  ampo .+= -t1,"Cdagup",b,"Cup",b+1
  ampo .+= -t1,"Cdagup",b+1,"Cup",b
  ampo .+= -t1,"Cdagdn",b,"Cdn",b+1
  ampo .+= -t1,"Cdagdn",b+1,"Cdn",b
  ampo .+= V1,"Ntot",b,"Ntot",b+1
end
for b=1:N-2
  ampo .+= -t2,"Cdagup",b,"Cup",b+2
  ampo .+= -t2,"Cdagup",b+2,"Cup",b
  ampo .+= -t2,"Cdagdn",b,"Cdn",b+2
  ampo .+= -t2,"Cdagdn",b+2,"Cdn",b
end
H = MPO(ampo,sites)

state = ["Emp" for n=1:N]
p = Ref(Npart)
for i=N:-1:1
  if p[] > i
    println("Doubly occupying site $i")
    state[i] = "UpDn"
    p[] -= 2
  elseif p[] > 0
    println("Singly occupying site $i")
    state[i] = (isodd(i) ? "Up" : "Dn")
    p[] -= 1
  end
end
# Initialize wavefunction to be bond
# dimension 10 random MPS with number
# of particles the same as `state`
psi0 = randomMPS(sites,state,10)

# Check total number of particles:
@show flux(psi0)

# Start DMRG calculation:
energy,psi = dmrg(H,psi0,sweeps)

upd = fill(0.0,N)
dnd = fill(0.0,N)
for j=1:N
  orthogonalize!(psi,j)
  psidag_j = dag(prime(psi[j], "Site"))
  upd[j] = scalar(psidag_j * op(sites, "Nup", j) * psi[j])
  dnd[j] = scalar(psidag_j * op(sites, "Ndn", j) * psi[j])
end

println("Up Density:")
for j=1:N
  println("$j $(upd[j])")
end
println()

println("Dn Density:")
for j=1:N
  println("$j $(dnd[j])")
end
println()

println("Total Density:")
for j=1:N
  println("$j $(upd[j]+dnd[j])")
end
println()

println("\nGround State Energy = $energy")


text:
include("Extensions.jl")

"""
### Route{T}
- path::String
- page::T \
A route is added to a ServerTemplate using either its constructor, or the
ServerTemplate.add(::Route) method. Each route calls either a particular
servable or function; the type of which denoted by T. The Route type is
    commonly constructed using the do syntax with the route(::Function, String)
    method.
##### example
```
# Constructors
route = Route("/", p(text = "hello"))

function example(c::Connection)
    write!(c, "hello")
end

route = Route("/", example)

# method
route = route("/") do c
    write!(c, "Hello world!")
    write!(c, p(text = "hello"))
    # we can also use extensions!
    c[:logger].log("hello world!")
end
```
------------------
##### field info
- **path::String**
The path, e.g. "/" at which to direct to the given component.
- **page::T** (::Function || T <: Component)
The servable to serve at this given route.
------------------
##### constructors
- Route(path::String, f::Function) **where**
- Route(path::String, s::Servable)
"""
mutable struct Route{T}
    path::String
    page::T
    function Route(path::String, f::Function)
        new{Function}(path, f)
    end
    function Route(path::String, s::Servable)
        new{typeof(s)}(path, s)
    end
end

"""
### ServerTemplate
- ip**::String**
- port**::Integer**
- routes**::Vector{Route}**
- extensions**::Dict**
- remove**::Function**
- add**::Function**
- start**::Function** \
The ServerTemplate is used to configure a server before
running. These are usually made and started inside of a main server file.
##### example

------------------
##### field info
- ip**::String**
- port**::Integer**
- routes**::Vector{Route}**
- extensions**::Dict**
- remove**::Function**
- add**::Function**
- start**::Function**
------------------
##### constructors
ServerTemplate(ip::String = "127.0.0.1", port::Int64 = 8001,
            routes::Dict = Vector{Route}()); extensions::Dict = Dict(:logger => Logger())
"""
mutable struct ServerTemplate
    ip::String
    port::Integer
    routes::Vector{Route}
    extensions::Dict
    remove::Function
    add::Function
    start::Function
    function ServerTemplate(ip::String = "127.0.0.1", port::Int64 = 8001,
        routes::Vector{Route} = Vector{Route}();
        extensions::Dict = Dict(:logger => Logger()))
        add, remove, start = serverfuncdefs(routes, ip, port, extensions)
        new(ip, port, routes, extensions, remove, add, start)::ServerTemplate
    end
end

"""
**Core**
### serverfuncdefs(::AbstractVector, ::String, ::Integer,
::Dict) -> (::Function, ::Function, ::Function)
------------------
This method is used internally by a constructor to generate the functions add,
start, and remove for the ServerTemplate.
#### example

"""
function serverfuncdefs(routes::AbstractVector, ip::String, port::Integer,
    extensions::Dict)
    add(r::Route{Function}) = push!(routes, r)
    add(r::Route{Servable}) = push!(routes, r)
    add(e::Any ...) = [push!(extensions, ext[1] => ext[2]) for ext in e]
    remove(i::Int64) = deleteat!(routes, i)
    start() = _start(routes, ip, port, extensions)
    return(add, remove, start)
end

"""
"""
mutable struct WebServer <: ToolipsServer
    host::String
    routes::Dict
    extensions::Dict
    server::Any
end

"""
**Core**
### _start(routes::AbstractVector, ip::String, port::Integer,
extensions::Dict) -> (::Sockets.HTTPServer)
------------------
This is an internal function for the ServerTemplate. This function is binded to
    the ServerTemplate.start field.
#### example

"""
function _start(routes::AbstractVector, ip::String, port::Integer,
     extensions::Dict)
    server = Sockets.listen(Sockets.InetAddr(parse(IPAddr, ip), port))
    if has_extension(extensions, Logger)
        extensions[Logger].log(1,
         "Toolips Server starting on port " * string(port))
    end
    routefunc, rdct, extensions = generate_router(routes, server, extensions)
    @async HTTP.listen(routefunc, ip, port, server = server)
    if has_extension(extensions, Logger)
        extensions[Logger].log(2,
         "Successfully started server on port " * string(port))
         extensions[Logger].log(1,
         "You may visit it now at http://" * string(ip) * ":" * string(port))
    end
    return(WebServer(ip, rdct, extensions, server))::WebServer
end

"""
**Core**
### generate_router(routes::AbstractVector, server::Any, extensions::Dict)
------------------
This method is used internally by the **_start** method. It returns a closure
function that both routes and calls functions.
#### example

"""
function generate_router(routes::AbstractVector, server, extensions::Dict)
    route_paths = Dict([route.path => route.page for route in routes])
    # Load Extensions
    ces::Dict = Dict()
    fes::Vector{ServerExtension} = Vector{ServerExtension}()
    for extension in extensions
        if typeof(extension[2].type) == Symbol
            if extension[2].type == :connection
                push!(ces, extension)
        elseif extension[2].type == :routing
                extension[2].f(route_paths, extensions)
            elseif extension[2].type == :func
                push!(fes, extension[2])
            end
        else
            if :connection in extension[2].type
                push!(ces, extension)
            end
            if :routing in extension[2].type
                extension[2].f(route_paths, extensions)
            end
            if :func in extension[2].type
                push!(fes, extension[2])
            end
        end
    end
    # Routing func

    routeserver::Function = function serve(http::HTTP.Stream)
        fullpath = http.message.target
        if contains(http.message.target, "?")
            fullpath = split(http.message.target, '?')[1]
        end
        c::Connection = Connection(route_paths, http, ces)
        if fullpath in keys(route_paths)
            if typeof(route_paths[fullpath]) <: Servable
                [extension.f(c) for extension in fes]
                route_paths[fullpath].f(c)
            else
                [extension.f(c) for extension in fes]
                route_paths[fullpath](c)
            end
        else
            if typeof(route_paths["404"]) <: Servable
                [extension.f(c) for extension in fes]
                route_paths[fullpath].f(c)
            else
                [extension.f(c) for extension in fes]
                route_paths["404"](c)
            end
        end
    end # serve()
    return(routeserver, route_paths, extensions)
end


text:
using DiffEqBase, OrdinaryDiffEq

type SimType{T} <: DEDataArray{T}
    x::Array{T,1}
    f1::T
end

f = function (t,u,du)
    du[1] = -0.5*u[1] + u.f1
    du[2] = -0.5*u[2]
end

const tstop1 = [5.]
const tstop2 = [8.]
const tstop = [5.;8.]

condition = function (t,u,integrator)
  t in tstop1
end

affect! = function (integrator)
  for c in user_cache(integrator)
    c.f1 = 1.5
  end
end

save_positions = (true,true)

cb = DiscreteCallback(condition, affect!, save_positions=save_positions)

condition2 = function (t,u,integrator)
  t in tstop2
end

affect2! = function (integrator)
  for c in user_cache(integrator)
    c.f1 = -1.5
  end
end

save_positions = (false,true)

cb2 = DiscreteCallback(condition2, affect2!, save_positions=save_positions)

cbs = CallbackSet(cb,cb2)

u0 = SimType{Float64}([10;10], 0.0)
prob = ODEProblem(f,u0,(0.0,10.0))
sol = solve(prob,Tsit5(),callback = cbs, tstops=tstop)

sol(1.5:0.5:2.5)

@test [sol[i].f1 for i in eachindex(sol)] == [zeros(9);1.5*ones(5);-1.5*ones(4)]


text:
###############################################################################
#
#  Maximal abelian subfield for fields function
#
###############################################################################

function check_abelian_extensions(class_fields::Vector{Tuple{ClassField{MapRayClassGrp, GrpAbFinGenMap}, Vector{GrpAbFinGenMap}}}, F::FieldsTower, IdExtension::GAP.GapObj)
  K = F.field
  autos = F.generators_of_automorphisms
  i = 1
  while codomain(F.subfields[i]) != K
    i += 1
  end
  #=
  E = GAP.Globals.SmallGroup(IdExtension)
  oE = GAP.gap_to_julia(Int, GAP.Globals.Size(E))
  #I need to compute the degree of the maximal abelian subextension over this subfield.
  deg_mas = Set{Int}()
  if degree(domain(F.subfields[i])) == 1
    idH = GAP.julia_to_gap([1, 1])
  else
    idH = IdGroup(automorphisms(domain(F.subfields[i])))
  end
  oH = GAP.gap_to_julia(Int, idH[1])
  order_sub = divexact(oE, oH)
  lN = GAP.Globals.NormalSubgroups(E)
  for s = 1:length(lN)
    if GAP.Globals.Size(lN[s]) != order_sub
      continue
    end
    if GAP.Globals.IdGroup(GAP.Globals.FactorGroup(E, lN[s])) == idH
      dS = GAP.gap_to_julia(Vector{Int}, GAP.Globals.AbelianInvariants(lN[s]))
      push!(deg_mas, prod(dS)
    end
  end
  return check_abelian_extensions(class_fields[indices], autos, F.subfields[i], deg_mas)
  =#
  return check_abelian_extensions(class_fields, autos, F.subfields[i])
end

function check_abelian_extensions(class_fields::Vector{Tuple{ClassField{MapRayClassGrp,GrpAbFinGenMap}, Vector{GrpAbFinGenMap}}}, autos::Vector{NfToNfMor}, emb_sub::NfToNfMor)#, deg_mas::Set{Int})
  @vprint :MaxAbExt 3 "Starting checking abelian extension\n"
  K = base_field(class_fields[1][1])
  d = degree(K)
  G = domain(class_fields[1][2][1])
  expo = G.snf[end]
  com, uncom = ppio(Int(expo), d)
  if com == 1
    return trues(length(class_fields))
  end
  #I extract the generators that restricted to the domain of emb_sub are the identity.
  #Notice that I can do this only because I know the way I am constructing the generators of the group.
  expG_arr = Int[]
  act_indices = Int[]
  p = 11
  d1 = discriminant(domain(emb_sub))
  d2 = discriminant(K)
  while iszero(mod(d1, p)) || iszero(mod(d2, p))
    p = next_prime(p)
  end
  R = GF(p, cached = false)
  Rx, x = PolynomialRing(R, "x", cached = false)
  fmod = Rx(K.pol)
  mp_pol = Rx(image_primitive_element(emb_sub))
  for i = 1:length(autos)
    pol = Rx(image_primitive_element(autos[i]))
    if mp_pol ==  Hecke.compose_mod(mp_pol, pol, fmod)
      push!(act_indices, i)
      #I compute the order of the automorphisms. I need the exponent of the relative extension!
      push!(expG_arr, _order(autos[i]))
    end
  end
  expG = lcm(expG_arr)
  expG1 = ppio(expG, com)[1]
  @vprint :MaxAbExt 3 "Context for ray class groups\n"

  OK = maximal_order(K)
  rcg_ctx = Hecke.rayclassgrp_ctx(OK, com*expG1)

  cfields = trues(length(class_fields))
  for i = 1:length(class_fields)
    @vprint :MaxAbExt 3 "Class Field $i\n"
    C, res_act = class_fields[i]
    res_act_new = Vector{GrpAbFinGenMap}(undef, length(act_indices))
    for i = 1:length(act_indices)
      res_act_new[i] = res_act[act_indices[i]]
    end
    cfields[i] = check_abelian_extension(C, res_act_new, emb_sub, rcg_ctx, expG)
  end
  return cfields

end

function check_abelian_extension(C::Hecke.ClassField, res_act::Vector{GrpAbFinGenMap}, emb_sub::NfToNfMor, rcg_ctx::Hecke.ctx_rayclassgrp, exponent_extension::Int)
  #I consider the action on every P-sylow and see if it is trivial
  G = codomain(C.quotientmap)
  expG = Int(exponent(G))
  fac = factor(rcg_ctx.n)
  prime_to_test = Int[]
  new_prime = false
  for (P, v) in fac
    # I check that the action of the P-Sylow has no fixed points.
    PS, mPS = psylow_subgroup(G, P, false)
    s, ms = snf(PS)
    act_sub = induce_action_on_subgroup(ms*mPS, res_act)
    if !isfixed_point_free(act_sub)
      push!(prime_to_test, P)
    end
  end
  if isempty(prime_to_test)
    return true
  end
  n = prod(prime_to_test)
  n1, m = ppio(Int(G.snf[end]), n)
  if m != 1
    # Now, I compute the maximal abelian subextension of the n-part of C
    Q, mQ = quo(G, n1, false)
    C1 = ray_class_field(C.rayclassgroupmap, C.quotientmap*mQ)
    #@vtime :MaxAbExt 1
    fl = _maximal_abelian_subfield(C1, emb_sub, rcg_ctx, exponent_extension)
  else
    #@vtime :MaxAbExt 1
    fl = _maximal_abelian_subfield(C, emb_sub, rcg_ctx, exponent_extension)
  end
  return fl

end

function _bound_exp_conductor_wild(O::NfOrd, n::Int, q::Int, bound::fmpz)
  d = degree(O)
  lp = prime_decomposition_type(O, q)
  f_times_r = divexact(d, lp[1][2])
  sq = fmpz(q)^f_times_r
  nbound = n+n*lp[1][2]*valuation(n,q)-div(fmpz(n), q^(valuation(n,q)))
  obound = flog(bound, sq)
  bound_max_ap = min(nbound, obound)  #bound on ap
  return div(q*bound_max_ap, n*(q-1)) #bound on the exponent in the conductor
end

function minimumd(D::Dict{NfOrdIdl, Int}, deg_ext::Int)
  primes_done = Int[]
  res = 1
  for (P, e) in D
    p = Int(minimum(P))
    if p in primes_done
      continue
    else
      push!(primes_done, p)
    end
    ram_index = P.splitting_type[1]
    s, t = divrem(e, ram_index)
    if iszero(t)
      d = min(s, valuation(deg_ext, p)+2)
      res *= p^d
    else
      d = min(s+1, valuation(deg_ext, p)+2)
      res *= p^d
    end
  end
  return res
end

function _maximal_abelian_subfield(A::Hecke.ClassField, mp::Hecke.NfToNfMor, ctx::Hecke.ctx_rayclassgrp, exp_extension::Int)
  K = base_field(A)
  k = domain(mp)
  ZK = maximal_order(K)
  zk = maximal_order(k)
  expected_order = Int(ppio(div(degree(K), degree(k)), ctx.n)[1])
  if gcd(expected_order, degree(A)) == 1
    return false
  end
  # disc(ZK/Q) = N(disc(ZK/zk)) * disc(zk)^deg
  # we need the disc ZK/k, well a conductor.
  d = abs(divexact(discriminant(ZK), discriminant(zk)^expected_order))
  mR1 = A.rayclassgroupmap
  expo = Int(exponent(codomain(A.quotientmap)))
  deg = expected_order
  #First, a suitable modulus for A over k
  #I take the discriminant K/k times the norm of the conductor A/K
  fm0 = Dict{NfOrdIdl, Int}()
  for (P, e) in mR1.fact_mod
    p = intersect_prime(mp, P)
    if !haskey(fm0, p)
      if !iscoprime(minimum(P, copy = false), expo)
        if e > 1
          fm0[p] = e
        end
      else
        fm0[p] = 1
      end
    end
  end
  lp = factor(ideal(zk, d))
  for (P, e) in lp
    if haskey(fm0, P)
      if !iscoprime(minimum(P, copy = false), deg*expo)
        fm0[P] = max(e, fm0[P])
      end
    else
      if !iscoprime(minimum(P, copy = false), deg*expo)
        fm0[P] = e
      else
        fm0[P] = 1
      end
    end
  end
  #Now, I extend this modulus to K
  fM0 = Dict{NfOrdIdl, Int}()
  for (p, v) in fm0
    lp = prime_decomposition(mp, p)
    if iscoprime(minimum(p, copy = false), expo*deg)
      for (P, e) in lp
        fM0[P] = 1
      end
    else
      for (P, e) in lp
        fM0[P] = e*v
      end
    end
  end
  ind = 0
  if isdefined(ctx, :computed)
    flinf = isempty(defining_modulus(mR1)[2])
    for i = 1:length(ctx.computed)
      idmr, ifmr, mRRR = ctx.computed[i]
      if flinf != ifmr
        continue
      end
      contained = true
      for (P, v) in fM0
        if !haskey(idmr, P) || idmr[P] < v
          contained = false
        end
      end
      if contained
        ind = i
        break
      end
    end
  end
  if iszero(ind)
    @vtime :MaxAbExt 1 R, mR = Hecke.ray_class_group_quo(ZK, fM0, defining_modulus(mR1)[2], ctx, check = false)
    if isdefined(ctx, :computed)
      push!(ctx.computed, (fM0, isempty(defining_modulus(mR1)[2]), mR))
    else
      ctx.computed = [(fM0, isempty(defining_modulus(mR1)[2]), mR)]
    end
  else
    mR = ctx.computed[ind][3]
    R = domain(mR)
  end
  if degree(zk) != 1
    if istotally_real(K) && isempty(defining_modulus(mR1)[2])
      inf_plc = InfPlc[]
    else
      inf_plc = real_places(k)
    end
    @vtime :MaxAbExt 1  r, mr = Hecke.ray_class_group(zk, fm0, inf_plc, n_quo = ctx.n)
  else
    rel_plc = true
    if istotally_real(K) && isempty(defining_modulus(mR1)[2])
      rel_plc = false
    end
    modulo = minimumd(fm0, expo * expected_order)
    wrp, trp = ppio(modulo, ctx.n)
    ftrp = factor(trp)
    for (pf, vpf) in ftrp
      if !iscoprime(pf-1, ctx.n)
        wrp *= pf
      end
    end
    @vtime :MaxAbExt 1  r, mr = Hecke.ray_class_groupQQ(zk, Int(wrp), rel_plc, ctx.n)
  end
  @vtime :MaxAbExt 1 lP, gS = Hecke.find_gens(mR, coprime_to = minimum(defining_modulus(mR1)[1]))
  listn = NfOrdIdl[norm(mp, x) for x in lP]
  # Create the map between R and r by taking norms
  preimgs = Vector{GrpAbFinGenElem}(undef, length(listn))
  for i = 1:length(preimgs)
    preimgs[i] = mr\listn[i]
  end
  proj = hom(gS, preimgs)
  #compute the norm group of A in R
  prms = Vector{GrpAbFinGenElem}(undef, length(lP))
  for i = 1:length(lP)
    if haskey(mR1.prime_ideal_preimage_cache, lP[i])
      f = mR1.prime_ideal_preimage_cache[lP[i]]
    else
      f = mR1\lP[i]
      mR1.prime_ideal_preimage_cache[lP[i]] = f
    end
    prms[i] = A.quotientmap(f)
  end
  proj1 = hom(gS, prms)
  S, mS = kernel(proj1, false)
  mS1 = mS*proj
  G, mG = Hecke.cokernel(mS1, false)
  return (order(G) == expected_order)::Bool

end


text:
const ALWB_URI = URI(scheme="http", host="www.bom.gov.au", path="/jsp/awra/thredds/fileServer/AWRACMS")

abstract type DataMode end

"""
    Values <: DataMode

Get as the regular measured values.
"""
struct Values <: DataMode end
"""
    Deciles <: DataMode

Get the dataset in relative deciles.
"""
struct Deciles <: DataMode end

struct ALWB{M<:DataMode,D<:Union{Day,Month,Year}} <: RasterDataSource end

layers(::Type{<:ALWB}) = (
    :rain_day, :s0_pct, :ss_pct, :sd_pct, :sm_pct, :qtot, :etot,
    :e0, :ma_wet, :pen_pet, :fao_pet, :asce_pet, :msl_wet, :dd
)

@doc """
    ALWB{Union{Deciles,Values},Union{Day,Month,Year}} <: RasterDataSource

Data from the Australian Landscape Water Balance (ALWB) data set.

See: [www.bom.gov.au/water/landscape](http://www.bom.gov.au/water/landscape)

Layers are available in daily, monthly and
annual resolutions, and as `Values` or relative `Deciles`.

The available layers are: `$(layers(ALWB))`.
""" ALWB

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/rain_day_2017.nc
# Precipiation = "rain_day"

# SoilMoisture_Upper = "s0_pct"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/ss_pct_2017.nc
# SoilMoisture_Lower = "ss_pct"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/sd_pct_2017.nc
# SoilMoisture_Deep = "sd_pct"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/sm_pct_2017.nc
# SoilMoisture_RootZone = "sm_pct"

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/qtot_2017.nc
# Runoff = "qtot"

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/etot_2017.nc
# Evapotrans_Actual = "etot"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/e0_2017.nc
# Evapotrans_Potential_Landscape = "e0"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/ma_wet_2017.nc
# Evapotrans_Potential_Areal = "ma_wet"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/pen_pet_2017.nc
# Evapotrans_Potential_SyntheticPan = "pen_pet"

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/fao_pet_2017.nc
# Evapotrans_RefCrop_Short = "fao_pet"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/asce_pet_2017.nc
# Evapotrans_RefCrop_Tall = "asce_pet"
# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/etot_2017.nc

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/msl_wet_2017.nc
# Evaporation_OpenWater = "msl_wet"

# http://www.bom.gov.au/jsp/awra/thredds/fileServer/AWRACMS/values/day/dd_2017.nc
# DeepDrainage = "dd"

"""
    getraster(T::Type{<:ALWB{Union{Deciles,Values},Union{Day,Month,Year}}}, layer; date)
    getraster(T::Type{<:ALWB{Union{Deciles,Values},Union{Day,Month,Year}}}, layer, date)

Download ALWB weather data, choosing layers from: `$(layers(ALWB))`.

Without a layer argument, all layers will be downloaded, and a tuple of path vectors returned.
If the data is already downloaded the path will be returned.
"""
getraster(T::Type{<:ALWB}, layer::Symbol; date) = getraster(T, layer, date)
function getraster(T::Type{<:ALWB{M,P}}, layer::Symbol, dates::Tuple) where {M,P}
    getraster(T, layer, _date_sequence(dates, P(1)))
end
function getraster(T::Type{<:ALWB}, layer::Symbol, date::Dates.TimeType)
    _check_layer(T, layer)
    mkpath(rasterpath(T))
    url = rasterurl(T, layer; date=date)
    path = rasterpath(T, layer; date=date)
    _maybe_download(url, path)
    path
end
function getraster(T::Type{<:ALWB}, layer::Symbol, dates::AbstractArray)
    getraster.(T, layer, dates)
end

rastername(T::Type{<:ALWB{M,P}}, layer; date) where {M,P} =
    string(layer, _pathsegment(P, date), ".nc")

rasterpath(::Type{ALWB}) = joinpath(rasterpath(), "ALWB")
rasterpath(::Type{ALWB{M,P}}) where {M,P} =
    joinpath(joinpath(rasterpath(), "ALWB"), map(_pathsegment, (M, P))...)
rasterpath(T::Type{<:ALWB}, layer; date=nothing) =
    joinpath(rasterpath(T), rastername(T, layer; date))

rasterurl(T::Type{<:ALWB{M,P}}, layer; date) where {M,P} =
    joinpath(ALWB_URI, _pathsegments(T)..., rastername(T, layer; date))

# Utilitiy methods

_pathsegments(::Type{ALWB{M,P}}) where {M,P} = _pathsegment(M), _pathsegment(P)
_pathsegment(::Type{Values}) = "values"
_pathsegment(::Type{Deciles}) = "deciles"
_pathsegment(::Type{Day}) = "day"
_pathsegment(::Type{Month}) = "month"
_pathsegment(::Type{Year}) = "year"
# Days are in whole-year files
_pathsegment(::Type{Day}, date) = "_" * string(year(date))
# Months and years are all in one file
_pathsegment(::Type{<:Union{Year,Month}}, date) = ""


text:
module SpringCollab2020RainbowSpinnerColorController

using ..Ahorn, Maple

@mapdef Entity "SpringCollab2020/RainbowSpinnerColorController" RainbowSpinnerColorController(x::Integer, y::Integer,
    colors::String="89E5AE,88E0E0,87A9DD,9887DB,D088E2", gradientSize::Number=280.0)

const placements = Ahorn.PlacementDict(
    "Rainbow Spinner Colour Controller (Spring Collab 2020)" => Ahorn.EntityPlacement(
        RainbowSpinnerColorController
    )
)

function Ahorn.selection(entity::RainbowSpinnerColorController)
    x, y = Ahorn.position(entity)

    return Ahorn.Rectangle(x - 12, y - 12, 24, 24)
end

Ahorn.render(ctx::Ahorn.Cairo.CairoContext, entity::RainbowSpinnerColorController, room::Maple.Room) = Ahorn.drawImage(ctx, Ahorn.Assets.northernLights, -12, -12)

end


max_stars_repo_name: pgrepds/Codegolf

max_issues_repo_name: pgrepds/Codegolf

max_forks_repo_name: pgrepds/Codegolf


text:
import Pkg;
import UUIDs;

offline_registry_name = Pkg.TOML.parsefile(
    joinpath(splitpath(@__DIR__)..., "Registry.toml")
    )["name"];
offline_registry_uuid_string = Pkg.TOML.parsefile(
    joinpath(splitpath(@__DIR__)..., "Registry.toml")
    )["uuid"];
offline_registry_uuid = UUIDs.UUID(offline_registry_uuid_string)
pushfirst!(Base.DEPOT_PATH,joinpath(splitpath(@__DIR__)..., "depot",),);
pushfirst!(
    Base.DEPOT_PATH,
    joinpath(
        homedir(),
        ".julia.isolated",
        offline_registry_name,
        offline_registry_uuid_string,
        ),
    );
filter!((x) -> !(Base.Filesystem.samefile(expanduser(x),expanduser("~/.julia")) || lowercase(strip(abspath(expanduser(x)))) == lowercase(strip(abspath(expanduser("~/.julia")))) || Base.Filesystem.samefile(x,"~/.julia") || lowercase(strip(abspath(x))) == lowercase(strip(abspath("~/.julia"))) || Base.Filesystem.samefile(expanduser(x),expanduser(joinpath(homedir(), ".julia"))) || lowercase(strip(abspath(expanduser(x)))) == lowercase(strip(abspath(expanduser(joinpath(homedir(), ".julia"))))) || Base.Filesystem.samefile(x,joinpath(homedir(), ".julia")) || lowercase(strip(abspath(x))) == lowercase(strip(abspath(joinpath(homedir(), ".julia"))))), Base.DEPOT_PATH);
unique!(Base.DEPOT_PATH);

try
    if !any(
            [
                (x.name == offline_registry_name &&
                    (x.uuid == offline_registry_uuid_string ||
                        x.uuid == offline_registry_uuid)) for
                            x in Pkg.Types.collect_registries()
                ]
            )
        Pkg.Registry.add(
            Pkg.RegistrySpec(path = joinpath(splitpath(@__DIR__)...),)
            );
    end
catch e
    @warn("ignoring exception: ", e,);
end

max_forks_repo_name: aluthge/OfflineRegistry