text:
pushopt(options::Dict,expr::Expr) = Base.Meta.isexpr(expr,:(=)) && (options[expr.args[1]] = expr.args[2])

const input_formats = @compat Dict{String, Any}(
        "noweb" => Dict{Symbol, Any}(
                    :codestart => r"^<<(.*?)>>=\s*$",
                    :codeend => r"^@\s*$"
                    ),
        "markdown" => Dict{Symbol, Any}(
                    :codestart => r"(?:^(?:`|~){3,}\s*(?:\{|\{\.|)julia(?:;|\s)(.*)\}\s*$)|(?:^(?:`|~){3,}\s*julia\s*$)",
                    :codeend => r"^`|~{3,}\s*$"
                    )
        )

function read_document(document, format)
  #doctext = readall(open(document))
  lines = split(bytestring(open(document) do io
                             mmap_array(Uint8,(filesize(document),),io)
                           end), "\n")

  codestart = input_formats[format][:codestart]
  codeend = input_formats[format][:codeend]
  state = "doc"

  docno = 1
  codeno = 1
  content = ""
  start_line = 0

  options = Dict()
  optionstring = ""
  parsed = Dict[]
  for lineno in 1:length(lines)
    line = lines[lineno]
    if (m = match(codestart, line)) != nothing && state=="doc"
      state = "code"
      if m.captures[1] == nothing
          optionstring = ""
      else
          optionstring=strip(m.captures[1])
      end
      #@show optionstring
      options = Dict{Symbol,Any}()
      if length(optionstring) > 0
          expr = parse(optionstring)
          Base.Meta.isexpr(expr,:(=)) && (options[expr.args[1]] = expr.args[2])
          Base.Meta.isexpr(expr,:toplevel) && map(pushopt,fill(options,length(expr.args)),expr.args)
      end
      haskey(options, :label) && (options[:name] = options[:label])
      haskey(options, :name) || (options[:name] = nothing)
      #@show options
      chunk = @compat Dict{Symbol,Any}(:type => "doc", :content => content,
                                       :number => docno,:start_line => start_line)
      docno += 1
      start_line = lineno
      push!(parsed, chunk)
      content = ""
      continue
    end
    if ismatch(codeend, line) && state=="code"
      chunk = @compat Dict{Symbol,Any}(:type => "code", :content => content,
                                       :number => codeno, :options => options,
                                       :optionstring => optionstring,
                                       :start_line => start_line)
      codeno+=1
      start_line = lineno
      content = ""
      state = "doc"
      push!(parsed, chunk)
      continue
    end

    content *= "\n" * line
  end

  #Remember the last chunk
  if content != ""
    chunk = @compat Dict{Symbol,Any}(:type => "doc", :content => content,
                                     :number =>  docno, :start_line => start_line)
    push!(parsed, chunk)
  end
  return parsed
end


text:
using BenchmarkTools
import LinearAlgebra, RecursiveFactorization

BenchmarkTools.DEFAULT_PARAMETERS.seconds = 0.5

luflop(m, n) = n^3÷3 - n÷3 + m*n^2
luflop(n) = luflop(n, n)

bas_mflops = Float64[]
rec_mflops = Float64[]
ns = 50:50:800
for n in ns
    A = rand(n, n)
    bt = @belapsed LinearAlgebra.lu!($(copy(A)))
    rt = @belapsed RecursiveFactorization.lu!($(copy(A)))
    push!(bas_mflops, luflop(n)/bt/1e9)
    push!(rec_mflops, luflop(n)/rt/1e9)
end

using Plots
plt = plot(ns, bas_mflops, legend=:bottomright, lab="OpenBLAS", title="LU Factorization Benchmark", marker=:auto, dpi=150)
plot!(plt, ns, rec_mflops, lab="RecursiveFactorization", marker=:auto)
xaxis!(plt, "size (N x N)")
yaxis!(plt, "GFLOPS")
savefig("lubench.png")
savefig("lubench.pdf")

ext: jl


text:
module IGC

using Dates
using TimeZones
using TimeZones: ZonedDateTime
import Base: parse, read, write, string, isempty

export IGCDocument, Abstract_IGC_record

include("formats.jl")
include("patterns.jl")
include("constants.jl")
include("exceptions.jl")
include("tlc.jl")
include("types.jl")
include("filename.jl")
include("reader.jl")
include("writer.jl")

end # module


text:
using Dates,StatsBase,DataFrames

include("constants.jl")

function match_dates(df::DataFrame,dates::Vector{Date};datecol::Symbol=:datetime)
    # TODO: check to make sure that dates range is within df[:,"datetime"] range
    # produces inconsistent results if not ... not sure why...
    mask = sum(map(x -> x.==df[:,datecol],dates),dims=1)[1]
    indices = findall(!iszero,mask)

    df[indices,:]

end

function upsample(img, dshape; method=mean)
    nshape = size(img)
    ydim,xdim = dshape

    window = [Int(round(v / dshape[i])) for (i,v) in enumerate(nshape)]
    wradius = Int.(round.(window / 2))

    outimg = Array{Float64,2}(undef,dshape)

    for i in 1:ydim, j in 1:xdim
        imin = (i-1) * window[1]
        imin = imin < 1 ? 1 : imin
        imax = imin + window[1]
        imax = imax > nshape[1] ? nshape[1] : imax
        jmin = (j-1) * window[2]
        jmin = jmin < 1 ? 1 : jmin
        jmax = jmin + window[2]
        jmax = jmax > nshape[2] ? nshape[2] : jmax

        outimg[i,j] = method(vec(img[imin:imax,jmin:jmax]))
    end

    return outimg
end

function dd_to_meters(lon::Float64, lat::Float64; scale::Float64=0.1)
    """
    Function to convert decimal degrees to meters based on the approximation
    given by: https://en.wikipedia.org/wiki/Geographic_coordinate_system
    Args:
    inPt (list or array): A Y,X point provided in geographic coordinates
    in that order.
    Keywords:
    scale (int): Resolution of the raster value to covert into meters,
    must be in decimal degrees.
    Returns:
    list: List of Y,X resolution values converted from meters to decimal degrees
    """

    radlat = deg2rad(lat) # convert degree latitude to radians

    ba::Float64 = 0.99664719 # constant of b/a

    ss = atan(ba*tan(radlat)) # calculate the reduced latitude

    # factor to convert meters to decimal degrees for X axis
    xfct = (pi / 180.0) * EARTHRADIUS * cos(ss)

    # factor to convert meters to decimal degrees for Y axis
    yfct = (111132.92-559.82 * cos(2*radlat)+1.175*cos(4*radlat)-
            0.0023*cos(6*radlat))

    # get meter resolution
    y_meters = scale * yfct
    x_meters = scale * xfct

    # return list of converted resolution values
    return hcat([x_meters], [y_meters])
end

function findnearest(lons,lats,x,y)
    xval,xidx = findmin(abs.(lons .- x))
    yval,yidx = findmin(abs.(lats .- y))

    return xidx,yidx
end

function meshgrid(x, y)
    X = [i for i in x, j in 1:length(y)]
    Y = [j for i in 1:length(x), j in y]
    return X, Y
end

function calc_area(lon,lat)
    dd_res = mean([abs(lon[1]-lon[2]),abs(lat[1]-lat[2])])
    lons,lats = meshgrid(lon,lat)
    sides = dd_to_meters.(lons,lats;scale=dd_res)

    return prod.(sides)
end


text:
using Cliquing
using Cliquing: greedycliquing!
using LinearAlgebra
using Memento
using Test

@testset "Cliquing.jl" begin
    @testset "Greedy Cliquing: Routine" begin
        @testset "AnyMaxClique" begin
            # AnyMaxClique
            a = Bool[1 1 0 0; 1 1 1 0; 0 1 1 1; 0 0 1 1]
            Cliquing.removediagonal!(a)
            out = Cliquing.anymaxclique(a)
            @test typeof(out) == GreedyClique
            @test Cliquing.mutuals(out) == BitVector([0, 0, 0, 0])
            @test Cliquing.vertices(out) == BitVector([0, 1, 1, 0])
        end

        @testset "NewClique" begin
            # NewClique
            connections = [true, false, true, false]
            v = [false, true, false, false]
            clique = GreedyClique(connections, v)
            @test Cliquing.mutuals(clique) == BitVector([1, 0, 1, 0])
            @test Cliquing.vertices(clique) == BitVector([0, 1, 0, 0])
        end

        @testset "CanMergeCliques" begin
            # CanMergeCliques
            c1 = GreedyClique(Bool[1, 0, 1, 0], Bool[0, 1, 0, 0])
            c2 = GreedyClique(Bool[0, 1, 0, 1], Bool[0, 0, 1, 0])
            out = Cliquing.mergeable(c1, c2)
            #@test size(out) == (4,)
            #@test out == [true, true, true, true]
            @test out == true
        end

        # UnionCliques
        @testset "UnionCliques" begin
            c1 = GreedyClique(Bool[1, 0, 1, 0], Bool[0, 1, 0, 0])
            c2 = GreedyClique(Bool[0, 1, 0, 1], Bool[0, 0, 1, 0])
            out = Cliquing.union(c1, c2)
            @test Cliquing.mutuals(out) == BitVector([0, 0, 0, 0])
            @test Cliquing.vertices(out) == BitVector([0, 1, 1, 0])
        end

        # GreedyCliquing
        @testset "GreedyCliqing" begin
            A = Bool[1 1 0 0; 1 1 1 0; 0 1 1 1; 0 0 1 1]
            cliques, singletons = greedycliquing(A, 2)
            @test A == A  # should not mutate input
            @test Cliquing.member(cliques[1]) == BitVector([0, 1, 1, 0])

            # Use the cliques [1,2,5] and [3,4] from https://en.wikipedia.org/wiki/Adjacency_matrix
            M = Bool[
                1 1 0 0 1 0
                1 0 1 0 1 0
                0 1 0 1 0 0
                0 0 1 0 1 1
                1 1 0 1 0 0
                0 0 0 1 0 0
            ]
            cliques, singletons = greedycliquing(M, 2)
            @test Cliquing.member(cliques[1]) == BitVector([1, 1, 0, 0, 1, 0])
            @test Cliquing.member(cliques[2]) == BitVector([0, 0, 1, 1, 0, 0])
            @test singletons == BitVector([0, 0, 0, 0, 0, 1])

            # Another test. Create a designed matrix where AnyMaxClique returns a clique
            # that does not contain the node with the largest number of neighbours
            a = Bool[
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
                0 0 0 0 1 0 0 0 0 0 0 0
            ]

            cliques, singletons = greedycliquing(a, 2)
            @test size(cliques) == (2,)
            @test Cliquing.member(cliques[2]) ==
                  BitVector([1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0])
            @test Cliquing.head(cliques[2]) ==
                  BitVector([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])
            @test Cliquing.member(cliques[1]) ==
                  BitVector([0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0])
            @test Cliquing.head(cliques[1]) ==
                  BitVector([0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0])
            @test singletons == BitVector([0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1])

            a = Bool[
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                1 1 1 1 0 0 0 0 0 0 0 0
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
                0 0 0 0 1 1 1 1 1 1 1 1
            ]

            cliques, singletons = greedycliquing(a, 2)
            @test size(cliques) == (2,)
            @test Cliquing.member(cliques[2]) ==
                  BitVector([1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0])
            @test Cliquing.head(cliques[2]) ==
                  BitVector([1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])
            @test Cliquing.member(cliques[1]) ==
                  BitVector([0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1])
            @test Cliquing.head(cliques[1]) ==
                  BitVector([0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0])
            @test singletons == BitVector([0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0])

            # Non-symmetric matrix input
            A = Bool[1 1 1 0; 1 1 1 0; 0 1 1 1; 0 0 0 1]
            @test_throws getlogger(@__MODULE__) ArgumentError greedycliquing(A, 2)

            @testset "Symmetric type" begin
                A = Symmetric(Bool[1 1 0 0; 1 1 1 0; 0 1 1 1; 0 0 1 1])
                cliques, singletons = greedycliquing(A, 2)
                @test A == A  # should not mutate input
                @test Cliquing.member(cliques[1]) == BitVector([0, 1, 1, 0])
                # Cannot mutate Symmetric matrix to non-symmetric matrix
                @test_throws MethodError greedycliquing!(A, 2)
            end
        end
    end
end


text:
import Base: >>>
export wrapbehavior,
       button,
       slider,
       checkbox,
       radio,
       radiogroup,
       togglebutton,
       textinput,
       progress,
       spinner,
       paper,
       datepicker

# A widget can be coerced into a behavior
# by calling `wrapbehavior` on it.
abstract Widget <: Behavior

subscribe(x::Signal, w::Widget) =
    subscribe(x, wrapbehavior(w))

# Ambiguity
intent(c::Union{Sampler, Collector}, tile::Widget) =
    OuterListener(c, tile)

intent(i::Intent, w::Widget) =
    intent(i, wrapbehavior(w))

default_intent(w::Widget) =
    default_intent(wrapbehavior(w))

## Button

@api button => (Button <: Widget) begin
    doc("A button.")
    arg(label::Tile, doc="The button label.")
    kwarg(
        raised::Bool=false,
        doc="If set to true, the button appears raised from the plane of the page."
    )
    kwarg(
        disabled::Bool=false,
        doc="If set to true, the button is disabled and will not be clickable."
    )
    kwarg(
        noink::Bool=false,
        doc="If set to true, disables the ripple effect when clicked."
    )
end

render(b::Button, state) =
    Elem("paper-button", render(b.label, state);
        attributes=@d(
            :raised => boolattr(b.raised),
            :noink => boolattr(b.noink),
            :disabled => boolattr(b.disabled)
        )
    )

wrapbehavior(b::Button) =
    clickable(b)

## Slider

@api slider => (Slider <: Widget) begin
    doc("""A slider. Use this to select values from within a continous range of
           numbers.""")
    arg(
        range::Range,
        doc="""The range specifying the minimum and maximum values that the slider
             can take."""
    )
    kwarg(
        value::Real=first(range),
        doc="The initial value of the slider. Defaults to the first value in the range."
    )
    kwarg(
        editable::Bool=true,
        doc="""If set to true, shows an editable text box with the current value
        of the slider."""
    )
    kwarg(
        pin::Bool=false,
        doc="""If set to true, shows a pin containing the current value as you
        drag the slider."""
    )
    kwarg(disabled::Bool=false, doc="If set to true, the slider is disabled.")
    kwarg(
        secondaryprogress::Real=0,
        doc="Highlight the slider bar from the beginning to this value."
    )
    kwarg(
        immediate::Bool=false,
        doc="""If set to true, also trigger while the user is dragging the slider.
        Due to the JavaScript implementation, this causes repeats in the subscribed
        values when the text box is edited. Use `droprepeats` from `Reactive.jl`
        to remove them."""
    )
end

wrapbehavior(s::Slider) =
    intent(ToType{eltype(s.range)}(),
           s.immediate ? hasstates(s, triggers=Dict("value"=>"change", "immediateValue"=>"immediate-value-change")) : hasstate(s))

render(s::Slider, state) =
    Elem("paper-slider", attributes=@d(
        :min => first(s.range),
        :max => last(s.range),
        :step => step(s.range),
        :value => s.value,
        :editable => boolattr(s.editable),
        :pin => boolattr(s.pin),
        :disabled => boolattr(s.disabled),
        :secondaryProgress => s.secondaryprogress,
        )
    )

## Checkbox
@api checkbox => (Checkbox <: Widget) begin
    doc("A checkbox.")
    typedarg(value::Bool=false, doc="State of the checkbox.")
    arg(label::Tile="", doc="The label.")
    kwarg(
        disabled::Bool=false,
        doc="If set to true, the checkbox will be disabled."
    )
end

wrapbehavior(c::Checkbox) =
    intent(ToType{Bool}(),
        hasstate(c, attr="checked", trigger="change"))

render(c::Checkbox, state) =
    Elem("paper-checkbox",
        render(c.label, state),
        attributes = @d(
            :checked=>boolattr(c.value),
            :disabled=>boolattr(c.disabled),
        )
    )

## Toggle Button
@api togglebutton => (ToggleButton <: Widget) begin
    doc("A toggle button.")
    arg(value::Bool=false, doc="State of the toggle button.")
    kwarg(
        disabled::Bool=false,
        doc="If set to true, the toggle button will be disabled."
    )
    kwarg(
        toggles::Bool=true,
        doc="If set to false, the user will not be able to change the state."
    )
end

wrapbehavior(c::ToggleButton) =
    intent(ToType{Bool}(),
        hasstate(c, attr="checked", trigger="change"))

render(c::ToggleButton, state) =
    Elem("paper-toggle-button",
        attributes = @d(
            :checked=>boolattr(c.value),
            :disabled=>boolattr(c.disabled),
            :toggles=>boolattr(c.disabled),
        )
    )

## Text input

@api textinput => (TextInput <: Widget) begin
    doc("A text input box.")
    arg(value::AbstractString="", doc="The current content.")
    kwarg(label::AbstractString="", doc="The label.")
    kwarg(error::AbstractString="", doc="Error to display if invalid input is entered.")
    kwarg(
        floatinglabel::Bool=true,
        doc="If set to true, the label floats above the input field when the
        input field is non-empty."
    )
    kwarg(
        multiline::Bool=false,
        doc="If set to true, input can contain new lines."
    )
    kwarg(rows::Int=0, doc="(Only in multiline mode). Number of rows of text.")
    kwarg(
        maxrows::Int=0,
        doc="""(Only in multiline mode). Maximum number of rows the input field
        will expand to. More lines will make the text input scrollable. (no limit if set to 0)""")
    kwarg(maxlength::Int=0, doc="Set the maximum length of input text.")
    kwarg(minlength::Int=0, doc="Set the minimum length of input text (not available in multiline).")
    kwarg(
        charcounter::Bool=false,
        doc="If set to true, a character count is displayed below the input field."
    )
    kwarg(
        pattern::AbstractString="",
        doc=md"""Pattern of allowed inputs. The pattern must match the entire value,
              not just some subset. The regular expression language is the same
              as [JavaScript's]
              (https://developer.mozilla.org/en/docs/Web/JavaScript/Guide/Regular_Expressions)
              (not available in multiline mode)
              ."""
    )
    #kwarg(autovalidate::Bool=true, doc="")
    kwarg(
        disabled::Bool=false,
        doc="If set to true, the text input will be disabled."
    )
end

wrapbehavior(t::TextInput, event="input") =
    hasstate(t, attr="value", trigger=event) |>
        intent(ToType{AbstractString}())

render(t::TextInput, state) = begin
    if t.multiline
        if length(t.pattern) > 0
            warn(
                "Multi-line text input does not support pattern validation")
        end

        elem = Elem("paper-textarea")

        if t.rows > 0
            elem &= @d(:attributes => @d(:rows => t.rows))
        end

        if t.maxrows > 0
            elem &= @d(:attributes => @d(:maxRows => t.maxrows))
        end
    else

        elem = Elem("paper-input")
        if t.pattern != ""
            elem &= @d(:attributes => @d("pattern" => t.pattern))
        end
    end

    if t.maxlength > 0
        elem &= @d(:attributes => @d(:maxlength => t.maxlength))
    end

    if t.minlength > 0
        elem &= @d(:attributes => @d(:minlength => t.minlength))
    end

    elem & @d(:attributes => @d(
                 "value" => t.value,
                 "label" => t.label,
                 "error-message" => t.error,
                 "auto-validate" => boolattr(true),
                 "disabled" => boolattr(t.disabled),
                 "char-counter" => boolattr(t.charcounter),
                 "no-label-float" => boolattr(!t.floatinglabel),
                 )
             )
end

## Radio buttons

@api radio => (RadioButton <: Tile) begin
    doc(md"""A radio button. Usually many radio buttons are grouped in a
    `radio group`.""")
    arg(
        name::AbstractString,
        doc=md"A name. The output of a `radiogroup` is the name of the selected radio button."
    )
    curry(label::Tile, doc="The label.")
    kwarg(
        value::Bool=false,
        doc="Is this radio button selected. When using in a radiogroup set `selected` field in the group instead."
    )
    kwarg(
        toggles::Bool=false,
        doc="If set to true, the radio button allows de-selection by clicking again."
    )
    kwarg(
        disabled::Bool=false,
        doc="If set to true, the radio button will be disabled."
    )
end

render(r::RadioButton, state) =
    Elem("paper-radio-button",
        render(r.label, state),
        attributes=@d(
            :name=>r.name,
            :checked=>boolattr(r.value),
            :toggles=>boolattr(r.toggles),
            :disabled=>boolattr(r.disabled)
        )
    )

@api radiogroup => (RadioGroup <: Widget) begin
    doc("""A group of radio buttons. At any time, only one radio button in a group
    can be selected.""")
    arg(radios::AbstractArray, doc="A vector of radio buttons.")
    kwarg(selected::AbstractString="", doc=md"Name of the currently selected `radiobutton`")
end

wrapradio(x::RadioButton) = x
wrapradio(x) = begin
    radio(string(x), x)
end

render(r::RadioGroup, state) =
    Elem("paper-radio-group",
        [render(wrapradio(b), state) for b in r.radios],
        attributes=@d(
            :selected=>r.selected,
        )
    )

wrapbehavior(r::RadioGroup) =
    hasstate(r, attr="selected", trigger="paper-radio-group-changed") |>
        intent(ToType{AbstractString}())

@api selector => (Selector <: Widget) begin
    arg(items::AbstractArray)
    kwarg(selected::Int=1)
end

## Spinner

@api spinner => (Spinner <: Tile) begin
    doc("A spinner. Usually used to denote something is loading or underway.")
    arg(active::Bool=true, doc="If set to false, the spinner will disappear.")
end

render(s::Spinner, state) =
    Elem("paper-spinner", attributes=@d(:active=>boolattr(s.active)))

## Progress bar

@api progress => (ProgressBar <: Tile) begin
    doc("A progress bar.")
    arg(value::Real, doc="Current primary progress.")
    kwarg(
        secondaryprogress::Real=0,
        doc="The secondary progress displayed in a lighter color."
    )
end

render(p::ProgressBar, state) =
    Elem("paper-progress";
        attributes = @d(
            "value"=>p.value,
            "secondary-progress"=>p.secondaryprogress,
        )
    )

@api paper => (PaperShadow <: Tile) begin
    doc("Raise a tile above the plane of the page and create a realistic shadow.")
    arg(elevation::Int, doc="The level to raise to. Valid values are Integers 1 to 5.")
    curry(tile::Tile, doc="The tile to be raised.")
    kwarg(
        animated::Bool=true,
        doc=md"If set to true, changes to `z` will be animated."
    )
end

render(p::PaperShadow, state) =
    Elem("paper-material", render(p.tile, state),
        attributes=@d(:elevation=>p.elevation, :animated=>p.animated))

# Date picker

# TODO: Migrate date picker to Polymer 1.0

if VERSION < v"0.4.0-dev"
    using Dates
end

@api dateselection => (DateSelection <: Behavior) begin
    arg(tile::Tile)
end
render(d::DateSelection, state) =
    render(d.tile, state) << Elem("date-selection")

immutable DateIntent <: Intent end
default_intent(::DateSelection) = DateIntent()

@api datepicker => (DatePicker <: Widget) begin
    doc("A date picker.")
    arg(date::Date=today(), doc=md"The date. Requires the `Dates` module on Julia v0.3")
    kwarg(
        range::Range{Date}=Date("1971-01-01"):Date("2100-12-31"),
        doc="The range of selectable dates."
    )
end
render(d::DatePicker, state) =
    Elem(
        "paper-date-picker-two",
        value=string(d.date),
        attributes=@d(:min=>string(first(d.range)), :max=>string(last(d.range)))
    )
wrapbehavior(p::DatePicker) = dateselection(p)

# TODO: Interpret as bounds error if date exceeds range
interpret(::DateIntent, d) = begin
    Date(d["year"], d["month"], d["day"])
end


text:
module VisionModel

using Flux

# Models
include("alexnet.jl")
include("vgg.jl")
include("resnet.jl")
include("googlenet.jl")

export  alexnet,
        vgg11, vgg11bn, vgg13, vgg13bn, vgg16, vgg16bn, vgg19, vgg19bn,
        ResNet18, ResNet34, ResNet50, ResNet101, ResNet152,
        googlenet

end


text:
module ReactiveMP

include("rocket.jl")
include("macrohelpers.jl")
include("helpers.jl")
include("math.jl")

include("constraints/prod/prod.jl")
include("constraints/form/form.jl")

include("message.jl")
include("marginal.jl")
include("distributions.jl")

include("constraints/prod/prod_analytical.jl")
include("constraints/prod/prod_generic.jl")
include("constraints/prod/prod_preserve_type.jl")
include("constraints/prod/prod_final.jl")

include("constraints/form/form_unspecified.jl")
include("constraints/form/form_point_mass.jl")
include("constraints/form/form_sample_list.jl")

as_marginal(message::Message)  = Marginal(getdata(message), is_clamped(message), is_initial(message))
as_message(marginal::Marginal) = Message(getdata(marginal), is_clamped(marginal), is_initial(marginal))

include("variable.jl")
include("pipeline.jl")

include("actors/score.jl")

include("algebra/helpers.jl")
include("algebra/cholinv.jl")
include("algebra/companion_matrix.jl")
include("algebra/permutation_matrix.jl")
include("algebra/correction.jl")

include("approximations.jl")
include("approximations/gausshermite.jl")
include("approximations/gausslaguerre.jl")
include("approximations/sphericalradial.jl")
include("approximations/laplace.jl")
include("approximations/importance.jl")
include("approximations/optimizers.jl")

include("distributions/pointmass.jl")
include("distributions/uniform.jl")
include("distributions/gamma_shape_rate.jl")
include("distributions/gamma.jl")
include("distributions/gamma_shape_likelihood.jl")
include("distributions/categorical.jl")
include("distributions/matrix_dirichlet.jl")
include("distributions/dirichlet.jl")
include("distributions/beta.jl")
include("distributions/bernoulli.jl")
include("distributions/normal_mean_variance.jl")
include("distributions/normal_mean_precision.jl")
include("distributions/normal_weighted_mean_precision.jl")
include("distributions/mv_normal_mean_covariance.jl")
include("distributions/mv_normal_mean_precision.jl")
include("distributions/mv_normal_weighted_mean_precision.jl")
include("distributions/normal.jl")
include("distributions/exp_linear_quadratic.jl")
include("distributions/wishart.jl")
include("distributions/contingency.jl")
include("distributions/function.jl")
include("distributions/sample_list.jl")

# Equality node is a special case and needs to be included before random variable implementation
include("nodes/equality.jl")

include("variables/random.jl")
include("variables/constant.jl")
include("variables/data.jl")

include("pipeline/async.jl")
include("pipeline/discontinue.jl")
include("pipeline/logger.jl")
include("pipeline/scheduled.jl")
include("pipeline/vague.jl")

include("rule.jl")

include("score.jl")
include("node.jl")

include("score/variable.jl")
include("score/node.jl")

# Stochastic nodes
include("nodes/uninformative.jl")
include("nodes/normal_mean_variance.jl")
include("nodes/normal_mean_precision.jl")
include("nodes/mv_normal_mean_covariance.jl")
include("nodes/mv_normal_mean_precision.jl")
include("nodes/mv_normal_weighted_mean_precision.jl")
include("nodes/gamma.jl")
include("nodes/gamma_shape_rate.jl")
include("nodes/beta.jl")
include("nodes/categorical.jl")
include("nodes/matrix_dirichlet.jl")
include("nodes/dirichlet.jl")
include("nodes/bernoulli.jl")
include("nodes/gcv.jl")
include("nodes/kernel_gcv.jl")
include("nodes/wishart.jl")
include("nodes/normal_mixture.jl")
include("nodes/gamma_mixture.jl")
include("nodes/dot_product.jl")
include("nodes/transition.jl")
include("nodes/autoregressive.jl")
include("nodes/bifm.jl")
include("nodes/bifm_helper.jl")
include("nodes/probit.jl")
include("nodes/flow/flow.jl")

# Deterministic nodes
include("nodes/addition.jl")
include("nodes/subtraction.jl")
include("nodes/multiplication.jl")

include("rules/prototypes.jl")

include("actors/prior.jl")

include("model.jl")

include("fixes.jl")

end


text:
module Rewrite

include("core.jl")
include("theory.jl")
include("interface.jl")

include("theory/free.jl")
include("theory/c.jl")
include("theory/ac.jl")

include("util.jl")

end # module


text:
###############################################################################
#
#   Cyclotomics
"""
    Cyclotomic(n, coeffs::AbstractVector)
Element of `n`-th cyclotomic field with coefficients stored as `coeffs`.

To access the internals of a cyclotomic use API functions:
 * `conductor` - the conductor of a cyclotomic, i.e. the `n` used currently for
 storage. This might not be the minimal embeding field of a cyclotomic.
 * `getindex`/`setindex!` - use `α[i]` to access the coefficient at `i`-th power
 of a cyclotomic (in a circular fashion)
 * `values`/`exponents` - paired iterators over _non zero_ coefficients/exponents
 corresponding to _non-zero_ coefficients
 * `normalform!` - bring a cyclotomic to its unique representation as given by
 Zumbroich basis (also available in non-modifying form).

Iteration over non-zero coefficients in `Cyclotomic` is provided by `iterate`
which produces pairs `(exp, coeff)` of exponent and corresponding coefficient.
"""
struct Cyclotomic{T,A<:AbstractVector{T}} <: Number
    n::Int
    coeffs::A
end

Cyclotomic(v::V) where {V<:AbstractVector} = Cyclotomic{eltype(v),V}(length(v), v)
Cyclotomic{T}(α::Cyclotomic) where {T} = Cyclotomic(conductor(α), convert.(T, coeffs(α)))
Cyclotomic{T,V}(α::Cyclotomic) where {T,V} =
    Cyclotomic{T,V}(conductor(α), convert.(T, coeffs(α)))

Cyclotomic{T,V}(a::R) where {T,V,R<:Real} = Cyclotomic{T,V}(1, T[a])

Cyclotomic(c::Complex{T}) where T = Cyclotomic{T, SparseVector{T, Int}}(c)
Cyclotomic{T,V}(c::C) where {T,V,C<:Complex} =
    Cyclotomic{T,V}(real(c)) + E(4) * Cyclotomic{T,V}(imag(c))

"""
    E(n[, i=1])
Return the `i`-th power of `n`-th root of unity with sparse vector as storage.
"""
function E(n, i = 1)
    k = Primes.totient(n)
    i = (0 <= i < n ? i : mod(i, n))
    coeffs = sparsevec([i + 1], [1], n)
    sizehint!(coeffs.nzind, k)
    sizehint!(coeffs.nzval, k)
    return Cyclotomic(n, coeffs)
end

####
#   Low-level interface

"""
    conductor(α::Cyclotomic)
Return the conductor, i.e. the degree of cyclotomic field `α` belongs to.
"""
conductor(α::Cyclotomic) = α.n

"""
    coeffs(α::Cyclotomic)
Return the coefficients of `α` as stored.
"""
coeffs(α::Cyclotomic) = α.coeffs

function _to_index(α::Cyclotomic, idx::Integer)
    # return mod(idx, conductor(α)) + 1
    0 <= idx < conductor(α) && return idx + 1
    conductor(α) <= idx && return (idx % conductor(α)) + 1
    return (idx % conductor(α)) + conductor(α) + 1
end

Base.@propagate_inbounds function Base.getindex(α::Cyclotomic, exp::Integer)
    return α.coeffs[_to_index(α, exp)]
end

Base.getindex(α::Cyclotomic, itr) = [α[i] for i in itr]

Base.@propagate_inbounds function Base.setindex!(α::Cyclotomic, val, exp::Integer)
    α.coeffs[_to_index(α, exp)] = val
    return val
end

struct ExpCoeffItr{C<:Cyclotomic}
    α::C
end

# general definitions for iteration
function Base.iterate(itr::ExpCoeffItr, state = 0)
    idx = findnext(!iszero, coeffs(itr.α), state + 1)
    idx === nothing && return nothing
    return (idx - 1, coeffs(itr.α)[idx]), idx
end

Base.IteratorSize(::Type{<:ExpCoeffItr}) = Base.HasLength()
Base.length(itr::ExpCoeffItr) = count(!iszero, coeffs(itr.α))
Base.eltype(::Type{<:ExpCoeffItr{<:Cyclotomic{T}}}) where {T} = Tuple{Int,T}

exps_coeffs(α::Cyclotomic) = ExpCoeffItr(α)
"""
    exponents(α::Cyclotomic)
Return an iterator over non-zero exponents of `α`, beginning at `0`-th one.
Matched iterator over coefficients is provided by @ref(values).
"""
exponents(α::Cyclotomic) = (e for (e, c) in exps_coeffs(α))

"""
    values(α::Cyclotomic)
Return an iterator over non-zero coefficients of `α`, beginning at `0`-th one.
Matched iterator over exponents is provided by @ref(exponents).
"""
Base.values(α::Cyclotomic) = (c for (e, c) in exps_coeffs(α))
Base.valtype(::Cyclotomic{T}) where {T} = T

Base.similar(α::Cyclotomic, T::Type = valtype(α)) = similar(α, T, conductor(α))
Base.similar(α::Cyclotomic, m::Integer) = similar(α, valtype(α), m)
Base.similar(α::Cyclotomic, T::Type, n::Integer) = Cyclotomic(similar(coeffs(α), T, n))

"""
    dense(α::Cyclotomic)
Return a copy of `α` with coefficients stored in dense `Vector`.
"""
dense(α::Cyclotomic{T}) where {T} = Cyclotomic{T,Vector{T}}(conductor(α), coeffs(α))

"""
    sparse(α::Cyclotomic)
Return a copy of `α` with coefficients stored in `SparseVector`.
"""
SparseArrays.sparse(α::Cyclotomic) = Cyclotomic(sparse(coeffs(α)))

function SparseArrays.droptol!(
    α::Cyclotomic{T,<:AbstractSparseVector},
    ε = eps(T) * 1e3,
) where {T}
    coeffs(α) .= droptol!(coeffs(α), ε)
    return α
end

function roundcoeffs!(
    α::Cyclotomic{<:AbstractFloat},
    r::Base.RoundingMode = Base.RoundNearest;
    kwargs...,
)
    for i in eachindex(coeffs(α))
        coeffs(α)[i] = round(coeffs(α)[i], r; kwargs...)
    end
    return α
end

function Base.Complex{T}(α::Cyclotomic) where {T<:AbstractFloat}
    z = zero(Complex{T})
    rα = reduced_embedding(α)
    n = conductor(rα)
    for (e, c) in exps_coeffs(rα)
        γ = 2 * T(π) * T(e) / n
        z += c * (cos(γ) + im * sin(γ))
    end
    return z
end

Base.complex(α::Cyclotomic{T}) where T = Complex{float(T)}(α)
Base.complex(α::Cyclotomic{T}) where T <: AbstractFloat = Complex{T}(α)

function Base.float(::Type{T}, α::Cyclotomic) where {T<:AbstractFloat}
    αre, αim = real.(Complex{T}.(reim(α)))
    if abs(αim/αre) <= sqrt(eps(T)) || αim < eps(T)
        return αre
    end
    throw(InexactError(:float, T, α))
end

Base.float(α::Cyclotomic) = float(Float64, α)
Base.Float64(α::Cyclotomic) = float(Float64, α)

function _isreal(α::Cyclotomic)
    rα = reduced_embedding(α)
    if !isreal(rα) || any(!iszero, (rα[i] for i in 1:conductor(rα)-1))
        return (false, rα)
    else
        return (true, rα)
    end
end

function Base.Int(α::Cyclotomic)
    flag, rα = _isreal(α)
    flag && return Int(rα[0])
    throw(InexactError(:Int, Int, α))
end

function Base.Rational{T}(α::Cyclotomic) where {T}
    flag, rα = _isreal(α)
    flag && return Rational{T}(rα[0])
    throw(InexactError(:Rational, Rational{T}, α))
end

Base.Rational(α::Cyclotomic{T}) where T<:Integer = Rational{T}(α)
Base.Rational(α::Cyclotomic{Rational{T}}) where T = Rational{T}(α)

Base.abs2(α::Cyclotomic) = α * conj(α)
Base.abs(α::Cyclotomic) = abs(complex(α))

Base.real(α::Cyclotomic) = (α + conj(α))/2
Base.imag(α::Cyclotomic) = -im*(α - conj(α))/2


text:
int_rules_6_3_12 = @theory begin

    #= ::Subsection::Closed:: =#
    #= 6.3.12*(d+e*x)^m*tanh(a+b*x+c*x^2)^n =#
    @apply_utils Antiderivative(tanh(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2) ^ ~(n'), ~x) => Integral(tanh(~a + ~b * ~x + ~c * (~x) ^ 2) ^ ~n, ~x) <-- FreeQ([~a, ~b, ~c, ~n], ~x)
    @apply_utils Antiderivative(coth(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2) ^ ~(n'), ~x) => Integral(coth(~a + ~b * ~x + ~c * (~x) ^ 2) ^ ~n, ~x) <-- FreeQ([~a, ~b, ~c, ~n], ~x)
    @apply_utils Antiderivative((~(d') + ~(e') * ~x) * tanh(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2), ~x) => (~e * log(cosh(~a + ~b * ~x + ~c * (~x) ^ 2))) / (2 * ~c) + ((2 * ~c * ~d - ~b * ~e) / (2 * ~c)) * Antiderivative(tanh(~a + ~b * ~x + ~c * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e], ~x)
    @apply_utils Antiderivative((~(d') + ~(e') * ~x) * coth(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2), ~x) => (~e * log(sinh(~a + ~b * ~x + ~c * (~x) ^ 2))) / (2 * ~c) + ((2 * ~c * ~d - ~b * ~e) / (2 * ~c)) * Antiderivative(coth(~a + ~b * ~x + ~c * (~x) ^ 2), ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e], ~x)
    #= Antiderivative((~x)^(~m)*tanh((~a')+(~b')*(~x)+(~c')*(~x)^2),~x) :=  x^(m-1)*log(cosh(a+b*x+c*x^2))/(2*c) - b/(2*c)*Antiderivative(x^(m-1)*tanh(a+b*x+c*x^2),x) - (m-1)/(2*c)*Antiderivative(x^(m-2)*log(cosh(a+b*x+c*x^2)),x) <--  FreeQ([a,b,c],x) && GtQ(m,1) =#
    #= Antiderivative((~x)^(~m)*coth((~a')+(~b')*(~x)+(~c')*(~x)^2),~x) :=  x^(m-1)*log(sinh(a+b*x+c*x^2))/(2*c) - b/(2*c)*Antiderivative(x^(m-1)*coth(a+b*x+c*x^2),x) - (m-1)/(2*c)*Antiderivative(x^(m-2)*log(sinh(a+b*x+c*x^2)),x) <--  FreeQ([a,b,c],x) && GtQ(m,1) =#
    @apply_utils Antiderivative((~(d') + ~(e') * ~x) ^ ~(m') * tanh(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2) ^ ~(n'), ~x) => Integral((~d + ~e * ~x) ^ ~m * tanh(~a + ~b * ~x + ~c * (~x) ^ 2) ^ ~n, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~m, ~n], ~x)
    @apply_utils Antiderivative((~(d') + ~(e') * ~x) ^ ~(m') * coth(~(a') + ~(b') * ~x + ~(c') * (~x) ^ 2) ^ ~(n'), ~x) => Integral((~d + ~e * ~x) ^ ~m * coth(~a + ~b * ~x + ~c * (~x) ^ 2) ^ ~n, ~x) <-- FreeQ([~a, ~b, ~c, ~d, ~e, ~m, ~n], ~x)
end

max_stars_repo_path: src/rules/6 Hyperbolic functions/6.3 Hyperbolic tangent/6.3.12 (d+e x)^m tanh(a+b x+c x^2)^n.jl


text:
function F_kendall(x, parms)
  X = x[1]
  α, μ, ν = parms

  birth = α * X
  death = μ * X
  immigration = ν

  return [birth, death, immigration]
end

function gillespie_kendall(;
  X = 5,
  α = 2.0, μ = 1.0, ν = 0.5)
  # initial conditions
  x0 = [X]
  # stoichiometries:
  #   birth, X --> X + X
  #   death, X --> 0
  #   immigration, 0 --> X
  nu = [1 -1 1]
  # parameter vector
  parms = [α, μ, ν]

  return x0, nu, parms
end

x0, nu, parms = gillespie_kendall()
F = F_kendall


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

using Test, Libdl
using LLVMLibUnwind_jll: llvmlibunwind_handle

@testset "LLVMLibUnwind_jll" begin
    if Sys.isapple()
        @test dlsym(llvmlibunwind_handle, :unw_getcontext; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_init_local; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_init_local_dwarf; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_step; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_get_reg; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_set_reg; throw_error=false) !== nothing
        @test dlsym(llvmlibunwind_handle, :unw_resume; throw_error=false) !== nothing
    end
end

max_stars_repo_name: jonas-schulze/julia

max_issues_repo_name: jonas-schulze/julia

max_forks_repo_name: jonas-schulze/julia


text:
#  Functions and methods common to all MixedModel types

"""
    feR(m::MixedModel)

Return th upper Cholesky factor for the fixed-effects parameters, as an `UpperTriangular`
`p × p` matrix.
"""
function feR(m::MixedModel)
    R = lmm(m).R
    kp1 = size(R, 1) - 1
    UpperTriangular(R[kp1, kp1])
end

"""
    lmm(m::MixedModel)

Extract the `LinearMixedModel` from a `MixedModel`.

If `m` is a `LinearMixedModel` return `m`. If `m` is a
`GeneralizedLinearMixedModel` return `m.LMM`.
"""
lmm(m::LinearMixedModel) = m
lmm(m::GeneralizedLinearMixedModel) = m.LMM

"""
    cond(m::MixedModel)

Returns the vector of the condition numbers of the blocks of `m.Λ`
"""
Base.cond(m::MixedModel) = [cond(λ)::Float64 for λ in lmm(m).Λ]

"""
    std{T}(m::MixedModel{T})

The estimated standard deviations of the variance components as a `Vector{Vector{T}}`.
"""
Base.std(m::MixedModel) = sdest(m) * push!([rowlengths(λ) for λ in lmm(m).Λ], [1.])

## methods for generics defined in StatsBase

function StatsBase.coeftable(m::MixedModel)
    fe = fixef(m)
    se = stderr(m)
    z = fe ./ se
    pvalue = ccdf(Chisq(1), abs2(z))
    CoefTable(hcat(fe, se, z, pvalue), ["Estimate", "Std.Error", "z value", "P(>|z|)"],
        coefnames(lmm(m).mf), 4)
end

"""
    describeblocks(io::IO, m::MixedModel)
    describeblocks(m::MixedModel)

Describe the types and sizes of the blocks in the upper triangle of `m.A` and `m.R`.
"""
function describeblocks(io::IO, m::MixedModel)
    lm = lmm(m)
    A, R = lm.A, lm.R
    for j in 1:size(A,2), i in 1:j
        println(io, i, ",", j, ": ", typeof(A[i,j]), " ", size(A[i,j]), " ", typeof(R[i,j]))
    end
    nothing
end
describeblocks(m::MixedModel) = describeblocks(Base.STDOUT, m)

"""
    fnames(m::MixedModel)

Returns the names of the grouping factors for the random-effects terms.
"""
function fnames(m::MixedModel)
    lm = lmm(m)
    [t.fnm for t in lm.wttrms[1:length(lm.Λ)]]
end

function getθ!(v::AbstractVector, m::MixedModel)
    Λ = lmm(m).Λ
    nl = map(nlower, Λ)
    if length(v) ≠ sum(nl)
        throw(DimensionMismatch("length(v) = $(length(v)) should be $(sum(nl))"))
    end
    offset = 0
    for i in eachindex(nl)
        nli = nl[i]
        getθ!(Compat.view(v, offset + (1 : nli)), Λ[i])
        offset += nli
    end
    v
end

function getθ(m::MixedModel)
    Λ = lmm(m).Λ
    getθ!(Array(eltype(Λ[1]), sum(A -> nlower(A), Λ)), m)
end

"""
    grplevels(m::MixedModel)

Returns the number of levels in the random-effects terms' grouping factors.
"""
function grplevels(m::MixedModel)
    lm = lmm(m)
    [length(lm.trms[i].f.pool) for i in eachindex(lm.Λ)]
end

"""
    ranef!{T}(v::Vector{Matrix{T}}, m::MixedModel{T}, uscale::Bool = false)

Overwrites `v` with the conditional modes of the random effects for `m`.

If `uscale` is `true` the random effects are on the spherical (i.e. `u`) scale, otherwise on the
original scale
"""
function ranef!(v::Vector, m::MixedModel, uscale)
    R, Λ = m.R, m.Λ
    k = length(Λ)        # number of random-effects terms
    for j in 1:k
        copy!(v[j], R[j, end])
    end
    rβ = R[k + 1, end]
    if !isempty(rβ)      #  in the pirls! function for GLMMs want to skip this
        β = vec(feR(m) \ rβ)
        kp1 = k + 1
        for j in 1:k     # subtract the fixed-effects contribution
            BLAS.gemv!('N', -1.0, R[j, kp1], β, 1.0, vec(v[j]))
        end
    end
    for j in k:-1:1
        Rjj = R[j, j]
        uj = vec(v[j])
        LinAlg.A_ldiv_B!(isa(Rjj, Diagonal) ? Rjj : UpperTriangular(Rjj), uj)
        for i in 1:j - 1
            ui = vec(v[i])
            ui -= R[i, j] * uj
        end
    end
    if !uscale
        for j in 1:k
            A_mul_B!(Λ[j], v[j])
        end
    end
    v
end

"""
    ranef{T}(m::MixedModel{T}, uscale=false)

Returns, as a `Vector{Matrix{T}}`, the conditional modes of the random effects in model `m`.

If `uscale` is `true` the random effects are on the spherical (i.e. `u`) scale, otherwise on the
original scale.
"""
function ranef(m::MixedModel, uscale=false)
    lm = lmm(m)
    Λ, trms = lm.Λ, lm.trms
    T = eltype(trms[end])
    v = Matrix{T}[]
    for i in eachindex(Λ)
        l = size(Λ[i], 1)
        k = size(trms[i], 2)
        push!(v, Array(T, (l, div(k, l))))
    end
    ranef!(v, lm, uscale)
end

"""
     vcov(m::MixedModel)

Returns the estimated covariance matrix of the fixed-effects estimator.
"""
function StatsBase.vcov(m::MixedModel)
    Rinv = inv(feR(m))
    varest(m) * (Rinv * Rinv')
end


text:
export CoordinateModel, CoordinateCreature

"""
# E.g.
minimizes function
    model = CoordinateModel(x -> abs(x[1] - 20.0), -200.0, 200.0)
    state = GAState(model, ngen=500, npop=6_000, elite_fraction=0.1,
                       crossover_rate=0.9, mutation_rate=0.9,
                       print_fitness_iter=1)
    ga(state)

    type T has to have properties
        y-x :: T
        0.25 .* (x+y) :: T
        randn(T) :: T
        z + 0.25*(y-x)*randn(T) :: T

"""
struct CoordinateModel{F,T} <: GAModel
    f::F
    xmin::T
    xmax::T
    xspan::T # xmax-xmin
    clamp::Bool
end
function CoordinateModel(f::F,xmin,xmax,clamp::Bool=true) where {F}
    xmin,xmax = promote(xmin,xmax)
    ET = eltype(xmin)
    N = length(xmin)
    xspan = xmax .- xmin
    # check that F(ymin), F(ymax) can be converted to Float64 without error
    z1 = Float64(f(xmin))
    z2 = Float64(f(xmax))
    # and that f(yspan), f(ymin), and f(ymax) has sane values maybe
    # z1!=Inf && z2!=Inf && !isnan(z1) && !isnan(z2) ||
    #    error("f(xmin) or f(xmax) objective function is either NaN or Inf")
    all(xspan .>= zero(ET)) || error("xmax[i] < xmin[i] for some i")
    CoordinateModel{F,typeof(xspan)}(f,xmin,xmax,xspan,clamp)
end

struct CoordinateCreature{T} <: GACreature
    value :: T
    objvalue :: Float64
end
CoordinateCreature(value::T, m::CoordinateModel{F,T}) where {F,T} =
    CoordinateCreature{T}(value, m.f(value))

fitness(x::CoordinateCreature{T}) where {T} = -x.objvalue

function randcreature(m::CoordinateModel{F,T}, aux, rng) where {F,T}
    xvalue = m.xmin .+ m.xspan .* rand(rng,length(m.xspan))
    CoordinateCreature(xvalue, m)
end

struct AverageCrossover end
function crossover(::AverageCrossover, z::CoordinateCreature{T},
                   x::CoordinateCreature{T}, y::CoordinateCreature{T},
                   m::CoordinateModel{F,T}, params, curgen::Integer,
                   aux, rng) where {F,T}
    z.value .= 0.5 .* (x.value .+ y.value)
    CoordinateCreature(z.value, m)
end

struct SinglePointCrossover end
function crossover(::SinglePointCrossover, z::CoordinateCreature{T},
                   x::CoordinateCreature{T}, y::CoordinateCreature{T},
                   m::CoordinateModel{F,T}, params, curgen::Integer,
                   aux, rng) where {F,T}
    N = length(x.value)
    i = rand(rng, 1:N)
    if rand(rng) < 0.5
        z.value[1:i] = x.value[1:i]
        z.value[i+1:end] = y.value[i+1:end]
    else
        z.value[1:i] = y.value[1:i]
        z.value[i+1:end] = x.value[i+1:end]
    end
    CoordinateCreature(z.value, m)
end

struct TwoPointCrossover end
function crossover(::TwoPointCrossover, z::CoordinateCreature{T},
                   x::CoordinateCreature{T}, y::CoordinateCreature{T},
                   m::CoordinateModel{F,T}, params, curgen::Integer,
                   aux, rng) where {F,T}
    N = length(x.value)
    i,j = rand(rng, 1:N, 2)
    i,j = i > j ? (j,i) : (i,j)
    if rand(rng) < 0.5
        z.value[:] = x.value
        z.value[i+1:j] = y.value[i+1:j]
    else
        z.value[:] = y.value
        z.value[i+1:j] = x.value[i+1:j]
    end
    CoordinateCreature(z.value, m)
end

function crossover(z::CoordinateCreature{T},
                   x::CoordinateCreature{T}, y::CoordinateCreature{T},
                   m::CoordinateModel{F,T}, params, curgen::Integer,
                   aux, rng) where {F,T}
    crossover(TwoPointCrossover(), z, x, y, m, nothing, curgen, aux, rng)
end

# Mutate over all dimensions
function mutatenormal(temp::Real, x::CoordinateCreature{T},
                      model::CoordinateModel{F,T}, rng) where {F,T}
    x.value .+= temp .* model.xspan .* randn(rng,length(x.value))
    model.clamp && (x.value .= clamp.(x.value, model.xmin, model.xmax))
    CoordinateCreature(x.value, model)
end

# Mutate through a single dimension
function mutatenormaldim(temp::Real, x::CoordinateCreature{T}, dim::Integer,
                         model::CoordinateModel{F,T}, rng) where {F,T}
    ET = eltype(T)
    x.value[dim] += temp * model.xspan[dim] * randn(rng,ET)
    model.clamp && (x.value[dim] = clamp(x.value[dim], model.xmin[dim], model.xmax[dim]))
    CoordinateCreature(x.value, model)
end

function mutate(x::CoordinateCreature{T}, model::CoordinateModel{F,T},
                params, curgen::Integer, aux, rng) where {F,T}
    if rand(rng) < params[:rate]
        if rand(rng) < get(params, :sa_rate, 0.0)
            sa(x,model,params[:k], params[:lambda],
               params[:maxiter], curgen, aux, rng)
        else
            N = length(x.value)
            mutatenormaldim(0.1, x, rand(1:N), model, rng)
        end
    else
        x
    end
end

export sa,satemp,saprob,mutatenormal
function sa(x::CoordinateCreature{T}, model::CoordinateModel{F,T},
            k::Real, lambda::Real, maxiter::Integer, curgen::Integer,
            aux, rng) where {F,T}
    N = length(x.value)
    y = x
    numnoups = 0 # number of consecutive
    #curgen = 0
    iter = curgen*maxiter + 0
    while true
        temp = satemp(iter, k, lambda) + numnoups/maxiter
        #temp = satemp(iter, k, lambda) + lambda * log(1+numnoups)
        dim = rand(rng,1:N)
        yvdim_old = y.value[dim]
        yov_old = y.objvalue
        fitness_old = fitness(y)
        y = mutatenormaldim(temp, y, dim, model, rng)
        diffe = fitness(y) - fitness_old
        if diffe >= 0
        elseif rand(rng) < saprob(diffe, temp)
        else
            y.value[dim] = yvdim_old
            y = CoordinateCreature(y.value, yov_old)
        end
        #println("temp: $temp newy: $(newy.objvalue) diffe: $diffe prob: $(saprob(diffe, temp))")
        #numnoups = ifelse(diffe > 0, 0, ifelse(diffe < 0, numnoups+1, numnoups))
        numnoups = ifelse(diffe > 0, 0, numnoups + 1)
        iter += 1
        iter > curgen*maxiter + maxiter && break
    end
    #println("pre: $(x.objvalue) post: $(y.objvalue)")
    y
end

satemp(iter::Integer, k::Real, lambda::Real) = k * exp(-lambda * iter)
# diff = newscore - oldscore
saprob(diff::Real, iter::Integer, k::Real, lambda::Real) =
    exp(diff / satemp(iter, k, lambda))
saprob(diff::Real, temp::Real) = exp(diff / temp)

selection(pop::Vector{<:CoordinateCreature{T}}, n::Integer, rng) where {T} =
    selection(TournamentSelection(2), pop, n, rng)

printfitness(curgen::Integer, x::CoordinateCreature{T}) where {T} =
    println("curgen: $curgen value: $(x.value) obj. value: $(x.objvalue)")


text:
using Test
using PlayingCards
using Logging
using TexasHoldem
using TexasHoldem: seat_number
TH = TexasHoldem

@testset "Table: constructors / observed cards" begin
    players = ntuple(i-> Player(Bot5050(), i), 2)
    deck = ordered_deck()
    shuffle!(deck)
    blinds = TH.Blinds(1,2)
    cards = TH.get_table_cards!(deck)
    table = TH.Table(;deck=deck, cards=cards, players=players)
    TH.deal!(table, blinds)

    table.state = PreFlop()
    @test TH.observed_cards(table) == ()
    table.state = Flop()
    @test TH.observed_cards(table) == table.cards[1:3]
    table.state = Turn()
    @test TH.observed_cards(table) == table.cards[1:4]
    table.state = River()
    @test TH.observed_cards(table) == table.cards
end

@testset "Table: Move button" begin
    # All players playing
    players = ntuple(i-> Player(Bot5050(), i), 3)
    table = Table(;players = players)
    @test TH.buttons(table.buttons) == (1, 2, 3, 1)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (2, 3, 1, 2)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (3, 1, 2, 3)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (1, 2, 3, 1)

    # Some players not playing
    players = (
        Player(Bot5050(), 1),
        Player(Bot5050(), 2; bank_roll=0),
        Player(Bot5050(), 3),
    )
    table = Table(;players = players)
    @test TH.buttons(table.buttons) == (1, 3, 1, 3)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (3, 1, 3, 1)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (1, 3, 1, 3)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (3, 1, 3, 1)

    # Some players not playing
    players = (
        Player(Bot5050(), 1),
        Player(Bot5050(), 2; bank_roll=0),
        Player(Bot5050(), 3),
    )
    table = Table(;players = players, dealer_id=2)
    @test TH.buttons(table.buttons) == (3, 1, 3, 1)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (1, 3, 1, 3)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (3, 1, 3, 1)
    move_buttons!(table)
    @test TH.buttons(table.buttons) == (1, 3, 1, 3)
end

@testset "Table: Circle index" begin
    @test TH.circle_index(5, 1) == 1
    @test TH.circle_index(5, 2) == 2
    @test TH.circle_index(5, 3) == 3
    @test TH.circle_index(5, 4) == 4
    @test TH.circle_index(5, 5) == 5
    @test TH.circle_index(5, 6) == 1
    @test TH.circle_index(5, 7) == 2
    @test TH.circle_index(5, 0) == 5
    @test TH.circle_index(5, -1) == 4
end

@testset "Table: Dealer iterator" begin
    @test TH.default_dealer_id() == 1

    # dealer_id = 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, Dealer(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [1, 2, 3, 4, 5]

    dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, Dealer(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [2, 3, 4, 5, 1]
end

@testset "Table: SmallBlind iterator" begin
    # dealer_id = 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, SmallBlind(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [2, 3, 4, 5, 1]

    # dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, SmallBlind(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [3, 4, 5, 1, 2]
end

@testset "Table: BigBlind iterator" begin
    # dealer_id = 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, BigBlind(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [3, 4, 5, 1, 2]

    # dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, BigBlind(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [4, 5, 1, 2, 3]
end

@testset "Table: FirstToAct iterator" begin
    # dealer_id = 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, FirstToAct(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [4, 5, 1, 2, 3]

    # dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, FirstToAct(), length(players)))
    sn = seat_number.(iter_players)
    @test sn == [5, 1, 2, 3, 4]
end

@testset "Table: iterate from player" begin
    # dealer_id = 1
    @test TH.default_dealer_id() == 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, players[1], length(players)))
    sn = seat_number.(iter_players)
    @test sn == [1, 2, 3, 4, 5]

    # dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, players[1], length(players)))
    sn = seat_number.(iter_players)
    @test sn == [1, 2, 3, 4, 5]

    # dealer_id = 1
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, players[2], length(players)))
    sn = seat_number.(iter_players)
    @test sn == [2, 3, 4, 5, 1]

    # dealer_id = 2
    players = ntuple(i-> Player(Bot5050(), i), 5)
    table = Table(;players = players, dealer_id = 2)
    TH.deal!(table, TH.blinds(table))
    iter_players = collect(TH.circle(table, players[2], length(players)))
    sn = seat_number.(iter_players)
    @test sn == [2, 3, 4, 5, 1]
end


text:
# LLR, "learning with linear rewards".
# Based on https://arxiv.org/abs/1011.4748.

mutable struct LLR <: Policy
  α::Float64
end

LLR() = LLR(1.0)

struct LLRDetails <: PolicyDetails
  solver_time::Float64
end

function choose_action(instance::CombinatorialInstance{T}, algo::LLR, state::State{T}; with_trace::Bool=false) where T
  if any(v == 0 for v in values(state.arm_counts))
    # Initialisation step: for each random variable (i.e. reward source), try arms, preferring those that have never been tested.
    weights = Dict(arm => (state.arm_counts[arm] == 0.0) ? 1.0 : 0.0 for arm in keys(state.arm_counts))
  else
    # Determine the weights for each arm and use them to solve the combinatorial problem.
    m = maximum_solution_length(instance)
    numerator = algo.α * sqrt((m + 1) * log(state.round))
    weights = Dict(arm => state.arm_average_reward[arm] + numerator * sqrt(1.0 / state.arm_counts[arm]) for arm in keys(state.arm_counts))
  end

  t0 = time_ns()
  sol = solve_linear(instance, weights)
  t1 = time_ns()

  if with_trace
    return sol, LLRDetails((t1 - t0) / 1_000_000_000)
  else
    return sol
  end
end


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

using Markdown
using InteractiveUtils

# ╔═╡ f84a00b0-3d55-11eb-1560-17af63da3f95
using Plots

# ╔═╡ e49da79c-3d55-11eb-2090-b1649653782c
begin
@userplot CirclePlot
@recipe function f(cp::CirclePlot)
    x, y, i = cp.args
    n = length(x)
    inds = circshift(1:n, 1 - i)
    linewidth --> range(0, 10, length = n)
    seriesalpha --> range(0, 1, length = n)
    aspect_ratio --> 1
    label --> false
    x[inds], y[inds]
end

n = 150
t = range(0, 2π, length = n)
x = sin.(t)
y = cos.(t)

anim = @animate for i ∈ 1:n
    circleplot(x, y, i)
end
gif(anim, "anim_fps15.gif", fps = 15)
end

# ╔═╡ Cell order:
# ╠═f84a00b0-3d55-11eb-1560-17af63da3f95
# ╠═e49da79c-3d55-11eb-2090-b1649653782c


text:
struct HeunEuler <: AbstractAdaptiveExplicitSolver
end

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

function n_intermediates(::HeunEuler)
    return 3 # to be checked
end

function n_temporaries(::HeunEuler)
    return 1
end

function stepper!(ode_solver::HeunEuler, ode_fcn, stepper_options, stepper_temp::StepperTemporaries)

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

    # views into the temporaries
    k1 = view(stepper_temp.der_x, 1:n_prob_dim, 1)
    k2 = view(stepper_temp.der_x, 1:n_prob_dim, 2)
    k3 = view(stepper_temp.der_x, 1:n_prob_dim, 3)

    x2 = view(stepper_temp.x_tmp, 1:n_prob_dim, 1)

    # intermediate steps
	t2 = t + 1 * dt
	@. x2 = x + 1 * dt * k1
	ode_fcn(k2, t2, x2)

    # new step
	stepper_temp.t_new    = t + dt
	@. stepper_temp.x_new = x + dt * (1/2 * k1 + 1/2 * k2)
	ode_fcn(k3, stepper_temp.t_new, stepper_temp.x_new)

    # error estimation
	@. stepper_temp.x_est = x + dt * 1 * k1
end

max_stars_repo_name: HomoModelicus/julia

max_issues_repo_name: HomoModelicus/julia

max_forks_repo_name: HomoModelicus/julia


text:
function invs(σx,σy; iσx,iσy,ms)
    iσx == 1 && iσy == 2 && return Invariants(ms,σ1=σx,σ2=σy)
    iσx == 2 && iσy == 1 && return Invariants(ms,σ1=σy,σ2=σx)
    iσx == 1 && iσy == 3 && return Invariants(ms,σ1=σx,σ3=σy)
    iσx == 3 && iσy == 1 && return Invariants(ms,σ1=σy,σ3=σx)
    iσx == 3 && iσy == 2 && return Invariants(ms,σ2=σy,σ3=σx)
    return Invariants(ms,σ2=σx,σ3=σy) # σx = 2 && σy = 3
    # error()
end

@recipe function f(ms::ThreeBodyMasses, intensity::Function)
    return (intensity, ms)
end

@recipe function f(intensity::Function, ms::ThreeBodyMasses; iσx = 1, iσy = 2, Ngrid = 100)
    #
    σxv = range(lims(iσx,ms)..., length=Ngrid)
    σyv = range(lims(iσy,ms)..., length=Ngrid)
    #
    #
    calv = [
        (Kibble(invs(σx,σy; iσx=iσx, iσy=iσy, ms=ms),ms^2) > 0 ?
        NaN :
        intensity(invs(σx,σy; iσx=iσx, iσy=iσy, ms=ms))) for σy in σyv, σx in σxv]

    seriestype := :heatmap
    colorbar --> false
    #
    σxv, σyv, calv
end

max_stars_repo_name: mmikhasenko/ThreeBodyDecay

max_issues_repo_name: mmikhasenko/ThreeBodyDecay

max_forks_repo_name: mmikhasenko/ThreeBodyDecay


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

# Collection of sources required to build MPCBuilder
name = "SymEngine"
version = v"0.4.0"
sources = [
    "https://github.com/symengine/symengine/releases/download/v$version/symengine-$version.tar.gz" =>
    "dd755901a9e2a49e53ba3bbe3f565f94265af05299e57a7b592186dd35916a1b",
]

# Bash recipe for building across all platforms

script = """cd symengine-$version
"""*raw"""
if [[ ${target} == x86_64-* ]]; then
    export CXXFLAGS=-march=core2
fi
cmake -DCMAKE_INSTALL_PREFIX=$prefix -DCMAKE_TOOLCHAIN_FILE=/opt/$target/$target.toolchain -DBUILD_TESTS=no -DBUILD_BENCHMARKS=no  -DBUILD_SHARED_LIBS=yes -DWITH_MPC=yes -DBUILD_FOR_DISTRIBUTION=yes -DWITH_SYMENGINE_THREAD_SAFE=yes .
make -j
make install
"""

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

# The products that we will ensure are always built
products(prefix) = [
    LibraryProduct(prefix, "libsymengine", :libsymengine)
]

# Dependencies that must be installed before this package can be built
dependencies = [
    "https://github.com/JuliaMath/GMPBuilder/releases/download/v6.1.2-2/build_GMP.v6.1.2.jl"
    "https://github.com/JuliaMath/MPFRBuilder/releases/download/v4.0.1-3/build_MPFR.v4.0.1.jl"
    "https://github.com/isuruf/MPCBuilder/releases/download/v1.1.0-2/build_MPC.v1.1.0.jl"
]

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

max_stars_repo_name: amcostant/SymEngineBuilder

max_issues_repo_name: amcostant/SymEngineBuilder

max_forks_repo_name: amcostant/SymEngineBuilder


text:
function dare(A, B, Q, R)
    n = size(A, 1)
    T = promote_type(eltype(A), eltype(B), eltype(Q), eltype(R))

    E = [
        Matrix{T}(I, n, n)  B/R*B'
        zeros(size(A))      A'
    ]

    F = [
        A   zeros(size(A))
        -Q  Matrix{T}(I, n, n)
    ]

    QZ = schur(F, E)
    QZ = ordschur(QZ, abs.(QZ.alpha./QZ.beta) .< 1)

    return QZ.Z[(n+1):end, 1:n]/QZ.Z[1:n, 1:n]
end

function care(A, B, Q, R)
    G = B/R*B'

    Z = [A   -G
         -Q  -A']

    S = schur(Z)
    S = ordschur(S, real(S.values) .< 0)
    U = S.Z

    (m, n) = size(U)
    U11 = U[1:div(m, 2), 1:div(n,2)]
    U21 = U[div(m,2)+1:m, 1:div(n,2)]

    return U21/U11
end

function lqr(A, B, Q, R)
    P = care(A, B, Q, R)
    K = R\B'*P
    return K
end

function dlqr(A, B, Q, R)
    P = dare(A, B, Q, R)
    K = (R + B'*P*B)\B'*P*A
    return K
end

function dlqr(A, B, Q, R, Δt)
    Q = Q*Δt
    R = R*Δt
    A = A*Δt+I
    B = B*Δt

    P = dare(A, B, Q, R)
    K = (R + B'*P*B)\B'*P*A
    return K
end


text:
include("Types.jl")

##=============判断是否可以行走=======================

# 判断是否可以向左走
function can_l(ma::maze, p)
  #colo = ma.arr[p...].color

  # 判断撞墙
  if p[1] - 1 < 1
    return false
  end

  # 判断撞到路径
  if ma.arr[(p + [-1, 0])...].color != 0
    return false
  end

  return true
end

function can_r(ma::maze, p)
  #colo = ma.arr[p...].color

  # 判断撞墙
  if p[1] + 1 > ma.m
    return false
  end

  # 判断撞到路径
  if ma.arr[(p + [1, 0])...].color != 0
    return false
  end

  return true
end

function can_u(ma::maze, p)
  #colo = ma.arr[p...].color

  # 判断撞墙
  if p[2] - 1 < 1
    return false
  end

  # 判断撞到路径
  if ma.arr[(p + [0, -1])...].color != 0
    return false
  end

  return true
end

function can_d(ma::maze, p)
  #colo = ma.arr[p...]

  # 判断撞墙
  if p[2] + 1 > ma.n
    return false
  end

  # 判断撞到自己
  if ma.arr[(p + [0, 1])...].color != 0
    return false
  end

  return true
end
#判断是否可以行走
function allCanMove(ma::maze, pa, pb)
  return (can_l(ma, pa) || can_r(ma, pa) || can_u(ma, pa) || can_d(ma, pa)) && (can_l(ma, pa) || can_r(ma, pa) || can_u(ma, pa) || can_d(ma, pa))
end

function can_move(ma::maze, p)
  return (can_l(ma, p) || can_r(ma, p) || can_u(ma, p) || can_d(ma, p))
end

##============尝试行走=========================

function nextp(ma::maze, p)
  t = false
  n = 0
  temp_p = p

  if !(can_l(ma, p) || can_r(ma, p) || can_u(ma, p) || can_d(ma, p))
    # 此点无法
    return p
  end

  while !t  && n < 1000
    n += 1
    k = rand([1, 2, 3, 4])
    if k == 1 && can_l(ma, p)#can_l
      t = true
      temp_p += [-1, 0]
      return temp_p
    end
    if k == 2 && can_r(ma, p)
      t = true
      temp_p += [1, 0]
      return temp_p
    end
    if k == 3 && can_u(ma, p)
      t = true
      temp_p += [0, -1]
      return temp_p
    end
    if k == 4 && can_d(ma, p)
      t = true
      temp_p += [0, 1]
      return temp_p
    end
  end
  throw("Cannot find any path here")
  return [-1, -1]
end

function finderr_makePathSuccess(ma::maze, pa, pb)
  if pa[1] < 1 || pa[1] > ma.m || pb[1] < 1 || pb[1] > ma.m
    throw("Break the bondary.")
  end

  if pa[2] < 1 || pa[2] > ma.n || pb[2] < 1 || pb[2] > ma.n
    throw("Break the bondary.")
  end
end

function is_near(pa, pb)
  d = sum(abs(pa - pb))
  return (d <= 1)
end

#回退
function walk_back(ma::maze, p, index)
  colo = index
  d = 20
  edx = 0

  if ma.p[index] > d
    edx = ma.p[index] - d
    if edx == 0
      throw("index = 0.")
    end
  else
    edx = 1
  end

  for i = edx + 1:ma.p[index]
    ma.arr[ma.path[index, i]...].color = 0
    ma.arr[ma.path[index, i]...].index = 0
  end
  ma.p[index] = edx
  if edx == 0
    throw("index = 0.")
  end

  return ma.path[index, edx]
end

#存入堆栈
function putin(ma::maze, t, point)
  ma.p[t] += 1
  ma.path[t, ma.p[t]] = point
  ma.arr[point...] = mycell(t, ma.p[t]) # 堆栈编号
end

# 判断周围是否有其他颜色的点
function get_near_one(ma, p)
  colo = ma.arr[p...].color
  co2 = 0
  around_move = reshape([1,0, -1,0, 0,1, 0,-1], 4, 2) # 四个方向

  for i = 1:4
    co2 = try
      ma.arr[p[1] + around_move[i,1], p[2] + around_move[i,2]].color
    catch
      -1
    end
    if co2 != colo && co2 != 0 && co2 != -1
      return ma.arr[p[1] + around_move[i,1], p[2] + around_move[i,2]].index
    end
  end

  # 未找到
  return -1
end

# 判断是否相碰，并以此修改堆栈
function get_near(ma, point_1, point_2)
  indexp = 0
  if (indexp = get_near_one(ma, point_1)) != -1
    # 清除其余路径
    for i = indexp + 1:ma.p[2]
      ma.arr[ma.path[2, i]...].color = 0
    end
    ma.p[2] = indexp  # 截断对面的路径
    return true
  elseif (indexp = get_near_one(ma, point_2)) != -1
    # 清除其余路径
    for i = indexp + 1:ma.p[1]
      ma.arr[ma.path[1, i]...].color = 0
    end
    ma.p[1] = indexp  # 截断路径
    return true
  end
  return false
end

# 尝试构建一次主干，成功则返回true,否则返回false
function makePathSuccess(ma::maze)
  pa = ma.st
  pb = ma.ed
  ma.p = [1, 1]
  ma.path[1,1] = pa
  ma.path[2,1] = pb
  ma.arr[pa...] = mycell(1, 1)
  ma.arr[pb...] = mycell(2, 1)
  n = 1
  #pa_old = pa

  #while allCanMove(ma, pa, pb) && n < 10000
  while (!get_near(ma, pa, pb)) && n < 100
    if can_move(ma, pa)
      pa = nextp(ma, pa)
      putin(ma, 1, pa) # 放入栈中
    else
      # 否则向回走
      pa = walk_back(ma, pa, 1)
      n += 1
    end

    if can_move(ma, pb)
      pb = nextp(ma, pb) # 各自颜色为1,2
      putin(ma, 2, pb)
    else
      pb = walk_back(ma, pb, 2)
      n += 1
    end

    # 染色
    #ma.arr[pa...].color = 1
    #ma.arr[pb...].color = 2
    #ma.arr[pa...].index = ma.p[1]
    #ma.arr[pb...].index = ma.p[2]

    # 防止越界
    finderr_makePathSuccess(ma, pa, pb)

  end

  #if n >= 10000
  #  throw("The path is too long")
  #end
  if get_near(ma, pa, pb)
    return true
  else
    #throw("Not near")
    return false
  end
end

function makePath(ma::maze, maxstep::Int)
  for i = 1:maxstep
    if makePathSuccess(ma)
      return 0
    else
      fill!(ma.arr, mycell(0, 0)) # 清除数据
      ma.p[1] = 1
      ma.p[2] = 1
    end
  end
  return -1
end

function reColorMaze(ma::maze)
  for i = 1:ma.m, j = 1:ma.n
    if ma.arr[i, j] != 0
      ma.arr[i, j] = 1
    end
  end
end

function makeMaze(m, n, st, ed, max_step::Int)
  temp = maze(fill(mycell(0, 0), m, n), m, n, st, ed, fill([0, 0], 2, m * n), [1, 1]) # 构建初始迷宫

  for i = 1:m, j = 1:n
    temp.arr[i, j].color = 0
    temp.arr[i, j].index = 0
  end

  temp.path[1, 1] = st
  temp.path[2, 1] = ed
  temp.arr[st...] = mycell(1,1)
  temp.arr[ed...] = mycell(2,1)

  #return temp

  #随机行走
  if makePath(temp, max_step) == -1
    throw("Cannot find any path in all trys.")
  end

  #reColorMaze(temp)

  return temp
end

#输出路径
function getPath(ma::maze)
  temp_path = Array(Int, sum(ma.p), 2)
  for i = 1:ma.p[1]
    temp_path[i, 1] = ma.path[1,i][1]
    temp_path[i, 2] = ma.path[1,i][2]
  end

  for j = 1:ma.p[2]
    temp_path[ma.p[1] + j, 1] = ma.path[2,ma.p[2] - j + 1][1]
    temp_path[ma.p[1] + j, 2] = ma.path[2,ma.p[2] - j + 1][2]
  end
  return temp_path
end

function getMap(ma::maze)
  m = ma.m
  n = ma.n

  temp_map = fill(0, m, n)
  for i = 1:m, j = 1:n
    temp_map[i, j] = ma.arr[i, j].color
  end
  return temp_map
end

max_stars_repo_name: RandomMaze/makeMaze

max_issues_repo_name: RandomMaze/makeMaze

max_forks_repo_name: RandomMaze/makeMaze


text:
## string search in the RegEx style

#=

API for pattern matching and replace operations redefining regular expression functions

The julia Regex type API consists of:

[1] *(r::Regex) in Base at regex.jl:677
[2] *(r1::Union{Regex, AbstractChar, AbstractString}, rs::Union{Regex, AbstractChar, AbstractString}...) in Base at regex.jl:656
[3] ==(a::Regex, b::Regex) in Base at regex.jl:614
[4] ^(r::Regex, i::Integer) in Base at regex.jl:729
[5] count(t::Union{Regex, AbstractString}, s::AbstractString; overlap) in Base at regex.jl:388
[6] eachmatch(re::Regex, str::AbstractString; overlap) in Base at regex.jl:609
[7] endswith(s::SubString, r::Regex) in Base at regex.jl:240
[8] endswith(s::AbstractString, r::Regex) in Base at regex.jl:235
[9] findall(t::Union{Regex, AbstractString}, s::AbstractString; overlap) in Base at regex.jl:361
[10] findfirst(r::Regex, s::AbstractString) in Base at regex.jl:327
[11] findnext(re::Regex, str::Union{String, SubString}, idx::Integer) in Base at regex.jl:299
[12] findnext(r::Regex, s::AbstractString, idx::Integer) in Base at regex.jl:324
[13] hash(r::Regex, h::UInt64) in Base at regex.jl:620
[14] match(re::Regex, str::Union{SubString{String}, String}, idx::Integer) in Base at regex.jl:274
[15] match(re::Regex, str::Union{SubString{String}, String}, idx::Integer, add_opts::UInt32) in Base at regex.jl:274
[16] match(r::Regex, s::AbstractString) in Base at regex.jl:294
[17] match(r::Regex, s::AbstractString, i::Integer) in Base at regex.jl:295
[18] occursin(r::Regex, s::SubString; offset) in Base at regex.jl:176
[19] occursin(r::Regex, s::AbstractString; offset) in Base at regex.jl:171
[20] show(io::IO, re::Regex) in Base at regex.jl:111
[21] startswith(s::SubString, r::Regex) in Base at regex.jl:208
[22] startswith(s::AbstractString, r::Regex) in Base at regex.jl:203

match(...) returns a MatchResult object.

Package Tokens implements similar methods for the following pattern types:

 1) search using Bloom algorithm, simple search and also extended search with multiple accepted bytes per byte position
 2) simple search using Boyer-Moore algorithm
 3) search for multiple patterns, each with syntax of 1)

Release 1.0: findnext, findfirst, startswith, endswith, occursin

replace:
 - simple search
 - bloom: jedes gefundene Zeichen kann für replace adressiert werden, Syntax [int:int] Zeichen von index (1 based) bis index.
   Negative int sind zulässig, -1 ist letztes byte -2 vorletztes ... [1:-1] ist also der gesamte Suchstring wie gefunden.

 - search syntax: ähnlich regex, muss auch uppercase ja/nein includieren

"""

"""
function Base.match(r::Regex, s::AbstractString, idx::Integer)
end

=#

"""
A Matcher supports string matching with precompiled patterns.

It consists of

 * a pattern definition in source code

 * a 'compiled' pattern structure (specific to an associated pattern matching algorithmn)

 * result of the last match attempt

A Matcher is not thread-safe and not reentrant.
It should be constructed and used in one thread only.
If matcher construction is expensive, it should support
the Base.copy function to create a copy with deep copies
of all fields which are mutable.

In the constructor, pattern definition is 'compiled' in a pattern structure
best suited for the associated pattern matching algorithm. In contrast to ordinary
string search, which does all search algorithm setup in every call, a Matcher
separates match initialization from the 'match kernel', which is some loop
over the positions in the string to search in.

The matching algorithm is implemented in the [`Base.match`](@ref) function, which
returns a match result object, and has a matcher as its first parameter.
All Matchers implemented in package Tokens return itself.

All search variants like findfirst, findnext, occursin, startswith, endswith and match
variants with several type varations on where to search are all mapped on the central
method

    Base.match(m::Matcher,ofs::UInt32, size::UInt64, s::String)

This is the only function (besides constructors) all matchers must implement.
Several signature variants are mapped with generic methods onto that central method.
This holds for Base.AbstractString, and AbstractToken and IOShared of this module.

Every concrete Matcher type must have the following properties:

* a field (or property) *base*::BToken which gives the byte sequence
  in which the search took place. Its token category is irrelevant
  (not used in matcher methods).

  * a field (or property) *matchsize*::UInt64 which is set to the number of
  bytes matched.

* a field (or property) *matchidx*::UInt32 which is set by a match call.
  If matching fails, it is set to 0. On match success,
  it is the index of the match within the buffer of *base*, and must fulfil
  offset(*base*) < *matchidx* <= offset(*base*)+usize32(*base*).
  If no match was found, *matchidx* must be equal to offset(*base*)
  its value is undefined.

Every concrete Matcher type which uses the generic search methods implemented
in package Tokens require additionally the following properties:

* a field (or property) *base*::BToken which gives the byte sequence
  in which the search took place. Its token category is irrelevant
  (not used in matcher methods)

"""
abstract type Matcher
end

# public Matcher API

# postponed [7] endswith(s::SubString, r::Matcher)
# postponed [8] endswith(s::AbstractString, r::Matcher)
# postponed [9] findall(t::Union{Matcher, AbstractString}, s::AbstractString; overlap)
# postponed [13] hash(r::Matcher, h::UInt64)
# postponed [14] match(re::Matcher, str::Union{SubString{String}, String}, idx::Integer)
# irrelevant[15] match(re::Matcher, str::Union{SubString{String}, String}, idx::Integer, add_opts::UInt32)
# postponed [21] startswith(s::SubString, r::Matcher)
# postponed [22] startswith(s::AbstractString, r::Matcher)

# setected API for release 1 of Tokens
#[10] findfirst(r::Matcher, s::AbstractString)
#[11] findnext(re::Matcher, str::Union{String, SubString}, idx::Integer)
#[12] findnext(r::Matcher, s::AbstractString, idx::Integer)
#[16] match(r::Matcher, s::AbstractString)
#[17] match(r::Matcher, s::AbstractString, i::Integer)
#[18] occursin(r::Matcher, s::SubString; offset)
#[19] occursin(r::Matcher, s::AbstractString; offset)
#[20] show(io::IO, re::Matcher)

Base.match(r::Matcher, s) = match(r,BToken(s))

Base.match(r::Matcher, s::BToken) = match(r,BToken(s))

# no. horrible because makes a copy ...
#Base.match(r::Matcher, s::AbstractString, start::Int) = match(r,string(s),start-firstindex(s)+1)

function Base.match(r::Matcher, s::Utf8String, start::Int)
    ss = substring(s)
    @boundscheck checkbounds(s,start)
    match(r,UInt32(start-1+offset(ss),usize(ss),ss.string))
    r.matchofs > 0 && (r.matchofs -= offset(ss))
    return r
end

function Base.match(r::Matcher, s::BToken, start::Int)
    @boundscheck checkbounds(s,start:last)
    match(r,start-1+offset(s),usize(s),s.buffer)
    r.matchofs > 0 && (r.matchofs -= offset(s))
    return r
end

Base.match(r::Matcher, s::AbstractToken, start::Int) = match(r,BToken(s),start)

Base.match(r::Matcher, io::IOShared) = match(r,io.readofs,usize(io),io.buffer)

"""
Simple string search structure using Bloom filters.
Algorithm is very similar to that used for string search in Julia Base.

Little initialization overhead, very fast, but no pattern flexibility.
"""
mutable struct BloomMatcher <: Matcher
  pattern :: String # pattern must match all bytes
  bloommask :: UInt64 # ORed byte hashes
  lastskip :: Int # bytes to skip if last byte matches but another byte does not
  bloomskip :: Int # bytes to skip if a byte-s hash is not in bloommask
  matchofs :: UInt32 #  offset (0-based) in string buffer to search of last match
  matchsize::UInt64 # 0 (not found) or number of bytes in match
  last :: UInt8 # last byte of pattern
  function BloomMatcher(pattern::AbstractString)
    matcher = new(string(pattern),0%UInt64,0,0,0,0x00)
    initialize!(matcher)
    return matcher
  end
end

"can be used to re-initialize after assigning a new pattern"
function initialize!(bm::BloomMatcher)
    p = bm.pattern
    n = ncodeunits(p)
    # n <= 1 && error("ExactMatcher requires a pattern size of at least 2")
    skip = n
    me.last = codeunit(p,n)
    bloom_mask = UInt64(_search_bloom_mask(tlast))
    bloom_bits = 1 # no. of bits set in bloom filter if <= 32
    bloom_skip = 1 # no. of bytes to skip if byte not in filter
    j = n
    while (j-=1)>=1  # j=n is already processed
        pj = codeunit(p,j)
        if pj == p.last && skip == n
            skip = n - j
        end
        if bloom_bits <= 32 # argument: is near max(p(bloom_skip)*bloom_skip)
            hash = bloomhash(pj)
            if hash&bloom_mask == 0
                if bloom_bits < 32
                    # put in bloom filter up to 32 bits
                    bloom_mask |= hash
                end
                bloom_bits += 1 # gets >32 to terminate bloom filter filling
            end
            if bloom_bits <= 32
                bloom_skip += 1
            end
        else
            if skip<n
                # we have finished adding hashes to bloom filter
                # and we have determined the skip distance if matching last byte
                # nothimg remains to be done in preprocessing - stop work.
                break
            end
        end
    end
    skip -= 1 # add 1 to skip in main loop
    bm.bloommask = bloom_mask
    bm.bloomskip = bloom_skip
    bm.lastskip = skip
end

"""
    match(matcher::ExactMatcher, s::String, sfirst::Int, slast::Int)

    s[sfirst:slast] is the string to be examined for pattern matches

    match returns matcher, all results are stored in matcher.

    boundscheck will ensure valid sfirst, slast are existing index values

"""
function Base.match(matcher::BloomMatcher,ofs::UInt32,size::UInt64, s::String)
    # TODO matchofs must be changed from index to offset, and matchsize must be set
    n = ncodeunits(matcher.pattern)
    @boundscheck checkrange(ofs,size,s)
    m = size%Int
    matcher.matchsize = 0 # assign default value "not found"
    if n<=1
        # special trivial case: switch to byte search or nothing to search
        if n==0
            m>0 && (matcher.matchofs = ofs+1)
        else
            # Base works like this: something(findnext(isequal(codeunit(e.last,1)), s, sfirst), 0)
            # let-s compare performance with an ordinary loop
            while sfirst <= slast
                @inbounds if matcher.last == codeunit(s,sfirst)
                    matcher.matchofs = sfirst
                end
                sfirst += 1
            end
            # result 0 is already assigned as 'default' - not necessary to assign 0
        end
    else
        w = m - n
        if w >= 0
            # w<0 implies: not found
            # ### START ### of pattern matching loop
            i = sfirst + n - 1 # position of last byte in matching candidate in search
            while i < m
                if codeunit(s,i) == matcher.last
                    # check candidate
                    j = 1
                    while j < n
                        if codeunit(s,i-n+j) != codeunit(matcher.pattern,j)
                            break
                        end
                        j += 1
                        # match found?
                        if j == n
                            matcher.matchofs = i-n+1
                            return matcher
                        end
                    end
                    # no match: skip and test bloom
                    i += matcher.lastskip
                    if i>=m # main loop finished?
                        break
                    end
                end
                i += 1
                if bloom_mask & _search_bloom_mask(codeunit(s,i)) == 0
                    i += matcher.bloomskip
                end
            end
            if i==m
                # test end match
                j = 1
                while j <= n
                    if codeunit(s,i-n+j) != codeunit(t,j)
                        break # not matchofs
                    end
                    j += 1
                end
                if j > n # means: loop without break , we have a match at the very end
                    matcher.matchofs = i-n+1
                end
           end
           # ###  END  ### of pattern matching loop
        end
    end
    return matcher
end

"""
Search for any byte from a given list.

"""
mutable struct AnyByteMatcher <: Matcher
    pattern :: Vector{UInt8} # pattern matches any of its bytes
    ismatch ::  Vector{UInt8} # 0 or index into pattern. index is 1+(byte value to test)
    matchofs :: UInt32 #  index (1-based) in string to search of last match or 0 (no match found)
    matchsize :: UInt64 # 0 or 1
    variant :: UInt32 # index into pattern for match
    function AnyByteMatcher(pattern::Vector{UInt8})
        ismatch =  Vector{UInt8}(0%UInt8,256)
        for i in 1::length(pattern)
            ismatch[pattern[i]+1] = i
        end
        new(pattern,ismatch,0,0)
    end
end

#function Base.match(matcher::AnyByteMatcher,s::String, sfirst::Int, slast::Int)
function Base.match(matcher::AnyByteMatcher,ofs::UInt32,size::UInt64, s::String)
    @boundscheck checkrange(ofs,size,s)
    @inbounds begin
    endofs = ofs + (size%UInt32)
        while ofs < endofs
            test = ismatch[byte(s,ofs)]
            if test> 0
                matcher.matchofs = ofs
                matcher.matchsize = 1
                matcher.variant = test
                return matcher
            end
            sfirst += 1
        end
        matcher.matchsize = 0
        return matcher
    end
end

"""
Restricted but fast RegEx matcher.

To achieve top search performance, it combines a Boyer-Moore alike
search which advances usually more than one position when no match
was found, with parallelized pattern tests using set operations.
It really boosts regex search performance for long, complicated
matches with lots of alternatives.

Matcher is constructed with an nonstandard string literal with
prefix *rr*. It is designed to be compatible with RegEx string constants,
implementing a RegEx subset: the quantifiers *, + and {..} are not supported.

String literal options are

    * i case insensitive matching (Ascii letters, only)

    * x white space (CR, LF, blank, tab) in the pattern definition
      is ignored, and everything from a '#' until end of line is
      skipped, allowing for comments.

    * l use linux style end-of-line encoding: a single LF.
      Affects the definitions of \\N and \$, see below

    * o use Apple OS style end-of-line encoding: a single CR.
      Affects the definitions of \\N and \$, see below

    * w use windows style end-of-line encoding: a CR LF sequence
      Affects the definitions of \\N and \$, see below

Most byte values in the pattern are treatet as exact match,
with the following exceptions:

    * \\x.. a byte value in hexadecimal notation, .. denotes two hex digits
      This is nonstandard for Regex implementations, but helpful if the
      matcher is working with binary patterns

    * \\d any ASCII digit (0..9).

    * \\D any byte which is NOT an ASCII digit (0..9)

    * \\w any ASCII word character (letter, digit, underscore).

    * \\W any byte which is NOT an ASCII word character (letter, digit, underscore)

    * \\w any ASCII word character (letter, digit, underscore).

    * \\W any byte which is NOT an ASCII word character (letter, digit, underscore)

    * \\s any white space byte: space, tab, newline, carriage return, vertical tab. With lexer: all continuations of raw category 16

    * \\S any byte which is NOT a white space byte: space, tab, newline, carriage return, vertical tab

    * \\r carriage return (0x13)

    * \\n line feed (0x10)

    * \\N neither \\r nor \\n

    * \\i switches for subsequent byte definitions to case unsensitive mode.
      Restriction: only ascii letters are covered, no further unicode letters.

    * \\I switches for subsequent byte definitions to case sensitive mode.
      This is the default, unless option 'i' is used.

    * * unsupported. Will cause a syntax error.
      In Regex: preceding byte/group may appear 0 or more times
      Workaround: try with a fixed number of repetitions like {0-5}

    * + unsupported. Will cause a syntax error.
      In Regex: preceding byte/group may appear 1 or more times
      Workaround: try with a fixed number of repetitions like {1-5}

    * ^ unsupported. Will cause a syntax error.
      In Regex: matches begin of data/line, depending on the option'm' set or unset.

    * \$ matches the end of line sequence, see options 'l', 'o', 'w'.
      In Regex: matches end of data/line, depending on the option 'm' set or unset.

    * . any byte
      In RegEx: same behavior, if option flag 's' is set.
      To match any byte which is not part of the end of line sequence, use \\N

    * \\ followed by any other byte: just that byte.
      Used to escape bytes which are otherwise interpreted as a command:
      \\ . * + [ ] ? ( ) { } ^ \$

    * [ ... ] matches one of the bytes defined in the list of byte definitions the brackets. All escapes above are allowed

    * [^...] any byte not in [...]

    * [.-.] matches bytes of a range. range is given by two byte definitions separated by '-'

    * [^.-.] any byte not in [.-.]

    * (...) a group definition. ... denotes a list of byte definitions.
      It may contain any escapes and also groups, recursively.
      The groups are numbered sequentially from left to right.
      After a match, all group values in the match can be requested as tokens.
      Example: the pattern "([a-z])([0-9])" matches "e5",
      group 1 value is then "e", group 2 value is "5".

      A group value can be the empty string if the group does not exist
      in the matching pattern alternative (but in other alternatives),
      or if the group is defined to be empty. Empty groups are used in
      pattern alternatives, see following chapters for ? and |.

    * \\1 .. \\9 Reference to a pattern group (first 9 groups).
      In a replacement string definition, it gets replaced by the group value.
      In a pattern, it will be replaced by the group definition
      (introducing an additional group), allowing re-use of lexer constructs.

    * ? preceding byte or group specification is optional in the match.
      Every ? is resolved by replacing the original pattern by two alternative
      patterns, one with the optional byte/group, one without it.
      In the latter, all groups are still there, but encoded as empty groups.
      Example: rr"((\\d)?\\d)" expands to two patterns, rr"((\\d)\\d)" and
      rr"(()\\d)". Group 2 is empty if only one digit was matched.

      n occurrences of '?' will result in 2^n patterns. Be aware of the technical
      limit in this implementation of up to 64 patterns.

    * {n} results in repeating the preceding byte or group n times.

    * {n-m} results in repeating the preceding byte or group
      a variable number of times, from n up to m. {0-1} is equivalent to ?.

      Processing is similar to ?: m-n+1 alternative patterns are created,
      one for each number from n to m. And if a group with inner groups
      is repeated, empty groups are added at the end of the repetitions,
      so that all variants have the same number of inner groups.

    * | separates alternative patterns. Can be used on top level
      or in groups.  Having m '|' in the pattern (or group) results in m+1
      alternative patterns. In each pattern, group index starts with the same value,
      there are no empty groups automatically added as in alternatives defined by '?'.
      If you want to distinguish groups coming from different alternatives, you must
      add sufficient many empty groups at the beginning of each alternative.

      m occurrences of '|' will result in m+1 patterns. This multiplies with the
      alternative patterns resulting from '?'. Be aware of the technical
      limit in this implementation of up to 64 alternative patterns after replacing
      all constructs defined by '?' and '|'.

    * {n} results in repeating the preceding byte or group n times.

    * {n-m} results in repeating the preceding byte or group
      a variable number of times, from n up to m. {0-1} is equivalent to ?.

A simple example for a date format: rs"(\\d?\\d)/\\1/\\1\\1"
It matches "1/2/21" and "01/02/2021", but also "99/0/389".
It expands to rr"(\\d?\\d)/(\\d?\\d)/(\\d?\\d)(\\d?\\d)",
showing that we have four groups. In a match, group 1 and 2 contain month
and day, group 3 and 4 concatenated give the year.

It is possible to define a matcher which accepts only valid day and month
numbers, and only two and four digit year numbers:
rr"((0?[123456789])|(1[012])):(\\2|([12]\\d)|(3[01])):(\\d\\d)?\\6".
Carefully counting the groups, we find the month in group 1 and also group 2,
the day in group 3 and also in group 4, the complete year in group 5,
century (if matched) in group 6 and the last two year digits in group 7.

There are still a few illegal day values accepted for months
with fewer days than 31, like "04/31/2021". Covering those
cases is possible, but requires different subpatterns for day
depending on the month, and correctly matching only valid dates
of the 29th of february is even dependent on the last two year digits.
"""
mutable struct RestrictedRegex <: Matcher
    pattern :: Vector{UInt8} # pattern matches any of its bytes
    ismatch ::  Vector{UInt8} # 0 or index into pattern. index is 1+(byte value to test)
    matchofs :: UInt32 #  index (1-based) in string to search of last match or 0 (no match found)
    variant :: UInt32 # index into pattern for match
    function RestrictedRegex(pattern,options)
      println("pattern size: ",ncodeunits(pattern))
      println("options: ",options)

        println("pattern: ",pattern)
        error("not yet implemented.")
    end
end

# TODO RegEx der Datum, Zeit, Int, Real erkennt. Datum: ca 10 Varianten. Zeit: ca. 5 Varianten (mit/ohne Sekunden).
# Zahlen bis 16 Ziffern, danach nix,. oder e/E => 48 Varianten. Parser danach muss selten int und immer real noch fertig parsen
# ggf auch 10 digits und 5 nach Dezimalpunkt, dann werden auch die meisten reals erwischt.

macro rr_str(arg,options)
    RestrictedRegex(arg,options)
end

AnyStringMatcher(pattern::AbstractVector{T} ) where T <: AbstractString = AnyStringMatcher(BTokenVector(pattern))

#=
function match(matcher::AnyStringMatcher,s::String, sfirst::Int, slast::Int )
    @boundscheck checkbounds(s,sfirst:slast)
    @inbounds begin
        while sfirst <= slast
            test = ismatch[codeunit(s,sfirst)]
            if test> 0
                matcher.matchofs = sfirst
                matcher.variant = test
                return matcher
            end
        end
        matcher.matchofs = 0
        return matcher
    end
end
=#


text:
export SimpleSparseBasis, SparseBasis, MaxBasis

# prototypes for space transforms and cutoffs
function transform end
function transform_d end
function transform_dd end
function inv_transform end

# ------------------------------------------------------------
#  Abstract polynomial degree business

"""
interface functions for `OneParticleBasis`
"""
function add_into_A! end

"""
interface functions for `OneParticleBasis`
"""
function add_into_A_dA! end

"""
`function degree` : compute some notion of degree of
the `arg` argument.
"""
function degree end

"""
`function scaling(b, p)`:

a scaling factor for a basis functions ϕ, which gives a rought estimate on
   the magnitude of ∇ᵖϕ e.g.,
```
ϕ = r^n Ylm
```
has scaling factor `n^p + l^p`, though sharper estimates are also possible.
"""
function scaling end

using LinearAlgebra: Diagonal

diagonal_regulariser(basis; diff = 0) = Diagonal(scaling(basis, diff))

"""
every scalar basis must implement this
"""
function rand_radial end


text:
using Documenter, AcqSynth

makedocs(
    sitename = "AcqSynth.jl",
    authors = "Jérémy Béjanin.",
    modules = [AcqSynth],
    linkcheck = true,
    clean = false,
    format = Documenter.HTML(prettyurls = false),
    pages = [
        "Home" => "index.md",
        "API" => [
            "api/helpers.md",
            "api/wrapper.md",
        ],
    ]
)

deploydocs(
    repo = "github.com/jebej/AcqSynth.jl.git",
    target = "build",
    branch = "gh-pages",
)


text:
# Compute the load in each edge based on the given routing and traffix matrix.
function compute_loads(routing::Routing, demands::Dict{Edge{Int}, Float64})
    rd = routing.data
    load = Dict{Edge{Int}, Float64}()
    for e in edges(rd)
        load[e] = 0.0
        for demand in keys(demands)
            if e in keys(routing.edge_flows[demand])
                load[e] += demands[demand] * routing.edge_flows[demand][e]
            end
        end
        load[e] /= capacity(rd, e)
    end
    return load
end

compute_loads(routing::RoutingSolution) =
    compute_loads(routing.routings[end])

compute_loads(routing::Routing) =
    compute_loads(routing, routing.data.traffic_matrix)

# Compute the max load based on the given routing and traffix matrix.
compute_max_load(routing::Routing, demands::Dict{Edge{Int}, Float64}) =
    maximum(values(compute_loads(routing, demands)))

compute_max_load(routing::RoutingSolution) =
    compute_max_load(routing.routings[end])

compute_max_load(routing::Routing) =
    compute_max_load(routing, routing.data.traffic_matrix)

# Compute the max load by the means of an optimisation problem.
compute_max_load(rd::RoutingData, demands::Dict{Edge{Int}, Float64}) =
    compute_max_load(_create_model(rd), rd, demands)

function compute_max_load(m::Model, rd::RoutingData, demands::Dict{Edge{Int}, Float64})
    rm = basic_routing_model_unitary(m, rd, rd.model_type.type)
    mu_capacity_constraints(rm, demands)
    @objective(m, Min, rm.mu)
    optimize!(m)
    return objective_value(m)
end