text:
function getfxpwidths(x::Real)
    xabs = abs(x)
    intwidth = Int(max(ceil(log2(floor(x))), 0))
    fracwidth = Int(max(ceil(-log2(xabs - floor(x))), 0))

    return (I = intwidth, F = fracwidth)
end

function getfixedpoint(x::Real, width)
    xabs = abs(x)
    intrep = Int(floor(xabs * 2^width[:F]))
    bitwidth = width[:I] + width[:F] + 1
    binstr = "$(x < 0 ? "-" : "")$bitwidth'b$(string(intrep; pad = bitwidth, base = 2))"

    return binstr
end

function constantreplacement!(m::Module, netlist::Netlist)
    maxintwidth = 1
    maxfracwidth = 0

    # find require widths
    for v in vertices(m.dfg)
        inputs = getinputs(m.dfg, v)
        outputs = getoutputs(m.dfg, v)
        op = Operation(gettype.(inputs), gettype.(outputs), getoperator(m.dfg, v))
        if allowconstreplace(gethandler(op))
            for input in inputs
                i = find(netlist, getname(input))
                if !isnothing(i) && isconstant(netlist[i])
                    width = getfxpwidths(parse(Float64, netlist[i].name))
                    maxintwidth = max(maxintwidth, width[:I])
                    maxfracwidth = max(maxfracwidth, width[:F])
                end
            end
        end
    end

    # replace all constants
    width = (I = maxintwidth, F = maxfracwidth)
    constreplacements = []
    for v in vertices(m.dfg)
        inputs = getinputs(m.dfg, v)
        outputs = getoutputs(m.dfg, v)
        op = Operation(gettype.(inputs), gettype.(outputs), getoperator(m.dfg, v))
        if allowconstreplace(gethandler(op))
            for (j, input) in enumerate(inputs)
                i = find(netlist, getname(input))
                if !isnothing(i) && isconstant(netlist[i])
                    conststr = getfixedpoint(parse(Float64, netlist[i].name), width)
                    push!(constreplacements, netlist[i].name => conststr)
                    inputs[j] = Variable(Symbol(conststr), input.type)
                end
            end
            set_prop!(m.dfg, v, :inputs, inputs)
        end
    end

    # update netlist
    for (orig, new) in unique(constreplacements)
        i = find(netlist, orig)
        update!(netlist, Net(name = new, class = :constant, size = (1, 1)))
    end

    return m
end


text:
@testset "CM types" begin
  Qx, x = QQ["x"]
  K, a = NumberField(x^4 + 8*x^2 + 4)
  cK = Hecke._complex_embeddings(K)
  @test length(cK) == 4
  E = cK[1]
  @test number_field(E) === K
  @test place(E) isa InfPlc
  @test overlaps(conj(E)(a), conj(E(a)))
  @test cK[1] == cK[1]
  @test cK[1] != cK[2]

  Ec = conj(E)
  fl, c = Hecke.iscm_field(K)
  @test c * E == Ec

  sprint(show, E)
  sprint(show, MIME"text/plain"(), E)

  k, b = NumberField(x^2 + 16)
  ktoK = hom(k, K, a^3 + 6*a)
  z = rand(k, -10:10)
  @test overlaps(restrict(E, ktoK)(z), E(ktoK(z)))

  ct = cm_types(k)
  cembd = Hecke._complex_embeddings(k)
  C = cm_type(k, cembd[1:1])
  @test number_field(C) === k
  @test Hecke.embeddings(C) == Hecke._complex_embeddings(k)[1:1]
  @test_throws ArgumentError cm_type(k, Hecke._complex_embeddings(k))
  @test_throws ArgumentError cm_type(number_field(x - 1)[1], Hecke.NumFieldEmbedding[])
  @test_throws ArgumentError cm_type(K, Hecke._complex_embeddings(K)[1:2]) # they are conjugated
  fl, c = Hecke.iscm_field(k)
  @test c * C == cm_type(k, cembd[2:2])
  @test id_hom(k) * C == C

  K, a = cyclotomic_field(7)
  ct = Hecke.cm_types(K)
  @test count(isprimitive, ct) == 6

  Qx, x = QQ["x"]
  f = x^4 + 104x^2 + 796 # primitive non-galois quartic CM field
  K, a = number_field(f)
  ct = Hecke.cm_types(K)
  KK = Hecke.reflex(ct[1]).field
  @test isisomorphic(KK, number_field(x^4 + 52*x^2 + 477)[1])[1]
  @test isisomorphic(Hecke.reflex(Hecke.reflex(ct[1])).field, K)[1] # since primitive
end

ext: jl


text:
using ConstituencyTrees, Test

@testset "ConstituencyTrees.jl" begin

    function testbrackets(t, str)
        buf = IOBuffer()
        Base.show(buf, t)
        result = String(take!(buf))
        ref = strip(replace(str, r"\s+"=>" "))
        @test result == "ConstituencyTree{String}" * ref
    end

    @testset "Trees" begin
        S = "(S (NP (DT the) (N cat)) (VP (V slept)))"
        tree = read_tree(S)
        @test tree.node == "S"

        @test collect(POS(tree)) == [("DT","the"),("N","cat"),("V","slept")]
        @test collect(Words(tree)) == ["the","cat","slept"]

        t2 = tree"(S (NP (DT the) (N cat)) (VP (V slept)))"
        @test tree == t2

        vals(x) = ConstituencyTrees.label.(collect(x))
        for (iter, nodes) in [(Leaves, "the cat slept"),
                              (PreOrderDFS, "S NP DT the N cat VP V slept"),
                              (PostOrderDFS, "the DT cat N NP slept V VP S")]
            @test vals(iter(tree)) == vals(iter(t2)) == split(nodes)
        end

        @test collect(Leaves(tree)) == collect(Leaves(t2)) == split("the cat slept")

        @test ConstituencyTrees.label(ConstituencyTree{String}("x")) == "x"
        @test ConstituencyTrees.isterminal(ConstituencyTree("DT", "the"))
        @test ConstituencyTrees.label(ConstituencyTree()) == nothing
        @test ConstituencyTrees.label("x") == "x"
        @test ConstituencyTrees.label(Tree(nothing)) == ""

        buf = IOBuffer()
        AbstractTrees.print_tree(buf, tree)
        @test String(take!(buf)) ==
            """
            S
            ├─ NP
            │  ├─ DT
            │  │  └─ "the"
            │  └─ N
            │     └─ "cat"
            └─ VP
               └─ V
                  └─ "slept"
            """
        @test ConstituencyTrees.height(tree) == 4
    end

    @testset "Sentiment" begin
        str = read(joinpath(@__DIR__, "data", "sentiment1.txt"), String)
        tree_str = read_tree(str)
        tree_int = read_tree(str, read_node = x -> parse(Int, x))
        @test collect(Words(tree_str)) == collect(Words(tree_int))
    end

    @testset "Brackets" begin
        S = "(S (NP (DT the) (N cat)) (VP (V slept)))"
        tree = read_tree(S)
    end

    @testset "Productions" begin
        tree = tree"(S (NP (DT the) (N cat)) (VP (V slept)))"
        @test productions(tree) == [
            ("S", ["NP","VP"]),
            ("NP", ["DT","N"]),
            ("DT", ["the"]),
            ("N", ["cat"]),
            ("VP", ["V"]),
            ("V", ["slept"])
        ]
    end

    @testset "Chomsky Normal Form" begin
        t = tree"(A B C D)"

        tree = tree"(A (B b) (C c) (D d))"

        testbrackets(chomsky_normal_form(tree, RightFactored()),
                     "(A (B b) (A|<C-D> (C c) (D d)))")

        ftree = tree"(A (B b) (C c) (D d) (E e) (F f))"
        testbrackets(chomsky_normal_form(ftree, LeftFactored()),
                     "(A (A|<B-C-D-E> (A|<B-C-D> (A|<B-C> (B b) (C c)) (D d)) (E e)) (F f))")

        testbrackets(chomsky_normal_form(ftree, RightFactored()),
                     "(A (B b) (A|<C-D-E-F> (C c) (A|<D-E-F> (D d) (A|<E-F> (E e) (F f)))))")

        pierre = tree"""
        (S
          (NP-SBJ
            (NP (NNP Pierre) (NNP Vinken))
            (, ,)
            (ADJP (NP (CD 61) (NNS years)) (JJ old))
            (, ,))
          (VP
            (MD will)
            (VP
              (VB join)
              (NP (DT the) (NN board))
              (PP-CLR
                (IN as)
                (NP
                  (DT a)
                  (JJ nonexecutive) (NN director)))
              (NP-TMP (NNP Nov.) (CD 29))))
          (. .))"""

        testbrackets(chomsky_normal_form(pierre, RightFactored()),
                     """
                     (S
                     (NP-SBJ
                         (NP (NNP Pierre) (NNP Vinken))
                         (NP-SBJ|<,-ADJP-,>
                         (, ,)
                         (NP-SBJ|<ADJP-,>
                             (ADJP (NP (CD 61) (NNS years)) (JJ old))
                             (, ,))))
                     (S|<VP-.>
                         (VP
                         (MD will)
                         (VP
                             (VB join)
                             (VP|<NP-PP-CLR-NP-TMP>
                             (NP (DT the) (NN board))
                             (VP|<PP-CLR-NP-TMP>
                                 (PP-CLR
                                 (IN as)
                                 (NP
                                     (DT a)
                                     (NP|<JJ-NN> (JJ nonexecutive) (NN director))))
                                 (NP-TMP (NNP Nov.) (CD 29))))))
                         (. .)))
                     """)

        testbrackets(chomsky_normal_form(pierre, LeftFactored()),
                     """
                     (S
                       (S|<NP-SBJ-VP>
                         (NP-SBJ
                           (NP-SBJ|<NP-,-ADJP>
                             (NP-SBJ|<NP-,> (NP (NNP Pierre) (NNP Vinken)) (, ,))
                             (ADJP (NP (CD 61) (NNS years)) (JJ old)))
                           (, ,))
                         (VP
                           (MD will)
                           (VP
                             (VP|<VB-NP-PP-CLR>
                               (VP|<VB-NP> (VB join) (NP (DT the) (NN board)))
                               (PP-CLR
                                 (IN as)
                                 (NP (NP|<DT-JJ> (DT a) (JJ nonexecutive)) (NN director))))
                             (NP-TMP (NNP Nov.) (CD 29)))))
                       (. .))
                     """)

        S = read(joinpath(@__DIR__, "data", "pierre.mrg"), String)
        tree = read_tree(S)
    end

    @testset "Collapse Unary" begin
        testbrackets(collapse_unary(tree"(A (A a))"), "(A (A a))")
        testbrackets(collapse_unary(tree"(A (A a))", collapse_root=true), "(A+A a)")
        testbrackets(collapse_unary(tree"((A (A a)))", collapse_root=false), "((A+A a))")

        pierre = tree"((S (NP-SBJ (NP (NNP Pierre) (NNP Vinken)) (, ,) (ADJP (NP (CD 61) (NNS years)) (JJ old)) (, ,)) (VP (MD will) (VP (VB join) (NP (DT the) (NN board)) (PP-CLR (IN as) (NP (DT a) (JJ nonexecutive) (NN director))) (NP-TMP (NNP Nov.) (CD 29)))) (. .)))"
        testbrackets(collapse_unary(pierre, collapse_root=false),
                     """((S
                            (NP-SBJ
                              (NP (NNP Pierre) (NNP Vinken))
                              (, ,)
                              (ADJP (NP (CD 61) (NNS years)) (JJ old))
                              (, ,))
                            (VP
                              (MD will)
                              (VP
                                (VB join)
                                (NP (DT the) (NN board))
                                (PP-CLR (IN as) (NP (DT a) (JJ nonexecutive) (NN director)))
                                (NP-TMP (NNP Nov.) (CD 29))))
                            (. .)))
                     """)
    end

    @testset "Treebanks" begin
        file = joinpath(@__DIR__, "data", "pierre.mrg")
        tb = Treebank(file)
        n = 0
        for tree in tb
            n += 1
        end
        @test n == 1
        @test length(collect(tb)) == 1

        tb2 = Treebank([file, file])
        @test length(collect(tb2)) == 2
        p1, p2 = Iterators.take(tb2, 2)
        @test p1 == p2

        tb = Treebank(joinpath(@__DIR__, "data"))
        @test length(collect(tb)) == 2
    end
end


text:
# Algorithm 2.1
# Bisection algorithm for finding the solutions to f(x) = 0

function bisect_approx(f, a, b, tolerance, maxIter)
	fa = f(a)
	for i ∈ 1:maxIter
		p = a + (b - a)/2
		fp = f(p)
		if fp == 0 || (b - a)/2 < tolerance
			return p
		end
		if sign(fa) == sign(fp)
			a = p
		else
			b = p
		end
	end
	throw(ArgumentError("Method failed"))
end

println(bisect_approx(sin, 1, 2π, 0.01, 100))


text:
# automatically generated by the FlatBuffers compiler, do not modify

MyGame.Example.eval(quote

struct AnyAmbiguousAliasesM1 end
AnyAmbiguousAliasesM1(args...; kwargs...) = Monster(args...; kwargs...)

struct AnyAmbiguousAliasesM2 end
AnyAmbiguousAliasesM2(args...; kwargs...) = Monster(args...; kwargs...)

struct AnyAmbiguousAliasesM3 end
AnyAmbiguousAliasesM3(args...; kwargs...) = Monster(args...; kwargs...)

FlatBuffers.@UNION(AnyAmbiguousAliases, (
    Nothing,
    AnyAmbiguousAliasesM1,
    AnyAmbiguousAliasesM2,
    AnyAmbiguousAliasesM3,
))

end)

max_stars_repo_name: rjkat/flatbuffers-julia

max_issues_repo_name: rjkat/flatbuffers-julia

max_forks_repo_name: rjkat/flatbuffers-julia


text:
function NamedData(arg0::JString, arg1::NamedData_StreamOpener)
    return NamedData((JString, NamedData_StreamOpener), arg0, arg1)
end

function get_name(obj::NamedData)
    return jcall(obj, "getName", JString, ())
end

function get_stream_opener(obj::NamedData)
    return jcall(obj, "getStreamOpener", NamedData_StreamOpener, ())
end


text:
function rrule(::Type{ITensor}, x1::AllowAlias, x2::TensorStorage, x3)
  y = ITensor(x1, x2, x3)
  function ITensor_pullback(ȳ)
    x̄1 = NoTangent()
    x̄2 = ȳ.tensor.storage
    x̄3 = NoTangent()
    return (NoTangent(), x̄1, x̄2, x̄3)
  end
  return y, ITensor_pullback
end

function rrule(::Type{ITensor}, x1::AllowAlias, x2::Tensor)
  y = ITensor(x1, x2)
  function ITensor_pullback(ȳ)
    x̄1 = NoTangent()
    x̄2 = Tensor(x1, ȳ)
    return (NoTangent(), x̄1, x̄2)
  end
  return y, ITensor_pullback
end

function rrule(::Type{ITensor}, x1::Tensor)
  y = ITensor(x1)
  function ITensor_pullback(ȳ)
    x̄1 = Tensor(ȳ)
    return (NoTangent(), x̄1)
  end
  return y, ITensor_pullback
end

function rrule(::typeof(itensor), x1::Tensor)
  y = itensor(x1)
  function itensor_pullback(ȳ)
    x̄1 = tensor(ȳ)
    return (NoTangent(), x̄1)
  end
  return y, itensor_pullback
end

function rrule(f::Type{<:Tensor}, x1::AllowAlias, x2::ITensor)
  y = f(x1, x2)
  function Tensor_pullback(ȳ)
    x̄1 = NoTangent()
    x̄2 = ITensor(x1, ȳ)
    return (NoTangent(), x̄1, x̄2)
  end
  return y, Tensor_pullback
end

function rrule(f::Type{<:Tensor}, x1::ITensor)
  y = f(x1)
  function Tensor_pullback(ȳ)
    x̄1 = ITensor(ȳ)
    return (NoTangent(), x̄1)
  end
  return y, Tensor_pullback
end

function rrule(::typeof(tensor), x1::ITensor)
  y = tensor(x1)
  function tensor_pullback(ȳ)
    x̄1 = ITensor(typeof(storage(x1))(ȳ.storage.data), inds(x1))
    return (NoTangent(), x̄1)
  end
  return y, tensor_pullback
end

# Special case for contracting a pair of ITensors
function ChainRulesCore.rrule(::typeof(*), x1::ITensor, x2::ITensor)
  y = x1 * x2
  function contract_pullback(ȳ)
    x̄1 = ȳ * dag(x2)
    x̄2 = dag(x1) * ȳ
    return (NoTangent(), x̄1, x̄2)
  end
  return y, contract_pullback
end

function ChainRulesCore.rrule(::typeof(*), x1::Number, x2::ITensor)
  y = x1 * x2
  function contract_pullback(ȳ)
    x̄1 = ȳ * dag(x2)
    x̄2 = dag(x1) * ȳ
    return (NoTangent(), x̄1[], x̄2)
  end
  return y, contract_pullback
end

function ChainRulesCore.rrule(::typeof(*), x1::ITensor, x2::Number)
  y = x1 * x2
  function contract_pullback(ȳ)
    x̄1 = ȳ * dag(x2)
    x̄2 = dag(x1) * ȳ
    return (NoTangent(), x̄1, x̄2[])
  end
  return y, contract_pullback
end

# TODO: use some contraction sequence optimization here
function ChainRulesCore.rrule(::typeof(*), x1::ITensor, x2::ITensor, xs::ITensor...)
  y = *(x1, x2, xs...)
  function contract_pullback(ȳ)
    tn = [x1, x2, xs...]
    N = length(tn)
    env_contracted = Vector{ITensor}(undef, N)
    for n in 1:length(tn)
      tn_left = tn[1:(n - 1)]
      # TODO: define contract([]) = ITensor(1.0)
      env_left = isempty(tn_left) ? ITensor(1.0) : contract(tn_left)
      tn_right = tn[reverse((n + 1):end)]
      env_right = isempty(tn_right) ? ITensor(1.0) : contract(tn_right)
      env_contracted[n] = dag(env_left) * ȳ * dag(env_right)
    end
    return (NoTangent(), env_contracted...)
  end
  return y, contract_pullback
end

function ChainRulesCore.rrule(::typeof(+), x1::ITensor, x2::ITensor)
  y = x1 + x2
  function add_pullback(ȳ)
    return (NoTangent(), ȳ, ȳ)
  end
  return y, add_pullback
end

function ChainRulesCore.rrule(::typeof(itensor), x::Array, a...)
  y = itensor(x, a...)
  function itensor_pullback(ȳ)
    uȳ = permute(unthunk(ȳ), a...)
    x̄ = reshape(array(uȳ), size(x))
    ā = broadcast_notangent(a)
    return (NoTangent(), x̄, ā...)
  end
  return y, itensor_pullback
end

function ChainRulesCore.rrule(::Type{ITensor}, x::Array{<:Number}, a...)
  y = ITensor(x, a...)
  function ITensor_pullback(ȳ)
    # TODO: define `Array(::ITensor)` directly
    uȳ = Array(unthunk(ȳ), a...)
    x̄ = reshape(uȳ, size(x))
    ā = broadcast_notangent(a)
    return (NoTangent(), x̄, ā...)
  end
  return y, ITensor_pullback
end

function ChainRulesCore.rrule(::Type{ITensor}, x::Number)
  y = ITensor(x)
  function ITensor_pullback(ȳ)
    x̄ = ȳ[]
    return (NoTangent(), x̄)
  end
  return y, ITensor_pullback
end

function ChainRulesCore.rrule(::typeof(dag), x)
  y = dag(x)
  function dag_pullback(ȳ)
    x̄ = dag(unthunk(ȳ))
    return (NoTangent(), x̄)
  end
  return y, dag_pullback
end

function ChainRulesCore.rrule(::typeof(permute), x::ITensor, a...)
  y = permute(x, a...)
  function permute_pullback(ȳ)
    x̄ = permute(unthunk(ȳ), inds(x))
    ā = broadcast_notangent(a)
    return (NoTangent(), x̄, ā...)
  end
  return y, permute_pullback
end

@non_differentiable combiner(::Indices)


text:
CURRENTDIR, CURRENTFILE = splitdir(@__FILE__)
ROOTDIR = splitdir(splitdir(CURRENTDIR)[1])[1]
OUTDIR = joinpath(ROOTDIR, "output", "two_planets")

# OUTDIR = "../../output/two_planets"

samplers = ["MALA", "AM", "SMMALA", "GAMC"]

base_speed = readcsv(joinpath(OUTDIR, "MALA", "summary.csv"))[end]

open(joinpath(OUTDIR, "summary.txt"), "w") do f
  for s in samplers
    out = readcsv(joinpath(OUTDIR, s, "summary.csv"))
    write(
      f,
      @sprintf(
        "%s & %.2f & %d & %d & %d & %d & %.2f & %.2f & %.2f\\\\\n",
        s,
        out[1],
        minimum(out[2:12]),
        mean(out[2:12]),
        median(out[2:12]),
        maximum(out[2:12]),
        out[13],
        out[end],
        out[end]/base_speed
      )
    )
  end
end


text:
function phylogeny_probabilities(
    replacement::Replacement, data::AbstractHLAData, tree::PhylogeneticTree
    )

    y = targets(replacement, data)
    ctree = deepcopy(tree)
    annotate!(ctree, data, replacement)

    p = state_probabilities(ctree)
    phylogeny_effect = [p[s]["1"] for s in string.(1:length(y))]
    phylogeny_effect = [min(max(0.01, x), 0.99) for x in phylogeny_effect]

    return phylogeny_effect
end

function phylogeny_probabilities(replacement::Replacement, data::AbstractHLAData)
    y = targets(replacement, data)
    tree = phylogenetic_tree(data)
    ctree = deepcopy(tree)
    annotate!(ctree, data, replacement)

    p = state_probabilities(ctree)
    phylogeny_effect = [p[s]["1"] for s in string.(1:length(y))]
    phylogeny_effect = [min(max(0.01, x), 0.99) for x in phylogeny_effect]

    return phylogeny_effect
end

function state_probabilities(tree::PhylogeneticTree)
    p = Dict{String, Dict{String, Float64}}()
    ctree = deepcopy(tree)

    states = ["0", "1"]

    observed_states = map(leaf -> get_property(ctree, leaf, :state), leaves(ctree))
    count_1 = count(observed_states .=== "1")
    count_0 = count(observed_states .=== "0")
    priors = Dict{String, Float64}(
        "0" => log(count_0 / (count_0 + count_1)),
        "1" => log(count_1 / (count_0 + count_1))
    )

    for leaf in leaves(ctree)
        leafname = get_property(ctree, leaf, :name)
        p[leafname] = Dict{String, Float64}()

        original_state = get_property(ctree, leaf, :state)
        loglikelihoods = Dict{String, Float64}()

        for state in states
            set_property!(ctree, leaf, :state, state)
            loglikelihoods[state] = L(ctree)
        end

        for state in states
            p[leafname][state] = exp((loglikelihoods[state]) -
                logsumexp((loglikelihoods[s]) for s in states))
        end

        set_property!(ctree, leaf, :state, original_state)
    end

    return p
end

function L(tree::PhylogeneticTree)
    fasta_file = write_fasta(tree)
    tree_file = write_newick(tree)

    output = read(`$(raxml_ng()) --evaluate --msa $fasta_file --tree $tree_file
        --model BIN --opt-model on --opt-branches off --threads 1 --nofiles`,
        String)

    ll = match(r"Final LogLikelihood: (.+)\n", output).captures[1]

    rm(fasta_file)
    rm(tree_file)

    return parse(Float64, ll)
end

function write_newick(tree::PhylogeneticTree)
    path, io = mktemp(cleanup = true)

    try
        write(io, newick_string(tree))
    finally
        close(io)
    end

    return path
end

function write_fasta(tree::PhylogeneticTree)
    path, io = mktemp(cleanup = true)

    try
        for (i, s) in enumerate(states(tree))
            write(io, ">$(i)\n")
            write(io, "$s\n")
        end
    finally
        close(io)
    end

    return path
end


text:
@testset "GeneralJONSWAP" begin
    θ₀ = [0.7, 0.9, 3.3, 5.0, 4.0]
    ω = 1.0

    @testset "Gradient" begin
        @test approx_gradient_uni(θ -> sdf(GeneralJONSWAP{1}(θ), ω), θ₀)      ≈ grad_sdf(GeneralJONSWAP{1}(θ₀), ω)
    end

    @testset "Error handling" begin
        @test_throws ArgumentError GeneralJONSWAP{1}(-0.7, 0.9, 3.3, 5.0, 4.0)
        @test_throws ArgumentError GeneralJONSWAP{1}(0.7, -0.9, 3.3, 5.0, 4.0)
        @test_throws ArgumentError GeneralJONSWAP{1}(0.7, 0.9, 0.8, 5.0, 4.0)
        @test_throws ArgumentError GeneralJONSWAP{1}(0.7, 0.9, 3.3, 0.9, 4.0)
        @test_throws ArgumentError GeneralJONSWAP{1}(0.7, 0.9, 3.3, 5.0, -1)
        @test_throws ArgumentError GeneralJONSWAP{1}(ones(4))
        @test_throws ArgumentError GeneralJONSWAP{1}(ones(6))
        @test_throws ErrorException GeneralJONSWAP(ones(5))
        @test_throws ErrorException GeneralJONSWAP(0.7, 0.9, 3.3, 5.0, 4.0)
    end
end


text:
#
# This file is a part of KronLinInv. License is MIT
# Copyright (c) 2019 Andrea Zunino
#

##======================================================================

@doc raw"""
    posteriormean_serial(klifac::KLIFactors,Gfwd::FwdOps,mprior::Array{Float64,1},
                  dobs::Array{Float64,1})

Computes the posterior mean model, purely serial version.

# Arguments
- `klifac`: a structure containing the required "factors" previously computed with
    the function `calcfactors()`. It includes
    * `U1, U2, U3` ``\mathbf{U}_1``, ``\mathbf{U}_2``, ``\mathbf{U}_3``  of  ``F_{\sf{A}}``
    * `diaginvlambda` ``F_{\sf{B}}``
    * `iUCmGtiCd1, iUCmGtiCd2, iUCmGtiCd3` ``\mathbf{U}_1^{-1}
      \mathbf{C}_{\rm{M}}^{\rm{x}}
      (\mathbf{G}^{\rm{x}})^{\sf{T}}(\mathbf{C}_{\rm{D}}^{\rm{x}})^{-1} ``,
      `` \mathbf{U}_2^{-1} \mathbf{C}_{\rm{M}}^{\rm{y}}
      (\mathbf{G}^{\rm{y}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{y}})^{-1} ``,
      `` \mathbf{U}_3^{-1} \mathbf{C}_{\rm{M}}^{\rm{z}}
      (\mathbf{G}^{\rm{z}})^{\sf{T}} (\mathbf{C}_{\rm{D}}^{\rm{z}})^{-1} ``
      of  ``F_{\sf{D}} ``
- `FwdOps`: a structure containing the three forward matrices
    * `G1, G2, G3` `` \mathbf{G} = \mathbf{G_1} \otimes \mathbf{G_2} \otimes \mathbf{G_3} ``
- `mprior`: prior model (vector)
- `dobs`:  observed data (vector)

# Returns
- The posterior mean model (vector)
"""
function posteriormean_serial(klifac::KLIFactors,Gfwd::FwdOps,mprior::Array{Float64,1},
                       dobs::Array{Float64,1})

    println("posteriormean_serial(): Serial version")

    ##--------------
    U1,U2,U3 = klifac.U1,klifac.U2,klifac.U3
    diaginvlambda = klifac.invlambda
    Z1,Z2,Z3 = klifac.iUCmGtiCd1,klifac.iUCmGtiCd2,klifac.iUCmGtiCd3
    G1,G2,G3 = Gfwd.G1,Gfwd.G2,Gfwd.G3

    ## sizes
    Ni = size(Z1,1)
    Nl = size(Z1,2)
    Nj = size(Z2,1)
    Nm = size(Z2,2)
    Nk = size(Z3,1)
    Nn = size(Z3,2)

    ## sizes
    # Nr12 = size(U1,2)*size(U2,2)*size(U3,2)
    # Nc1  = size(Z1,1)*size(Z2,1)*size(Z3,1)
    Na = size(mprior,1)
    Nb = size(dobs,1)

    ##-------------
    av = collect(1:Na)
    bv = collect(1:Nb)
    ## vectors containing all possible indices for
    ##    row calculations of Kron prod AxBxC
    iv = div.( (av.-1), (Nk.*Nj) ) .+1
    jv = div.( (av.-1 .-(iv.-1).*Nk.*Nj), Nk ) .+1
    kv = av.-(jv.-1).*Nk.-(iv.-1).*Nk.*Nj
    ## vectors containing all possible indices for
    ##    column calculations of Kron prod AxBxC
    lv =  div.( (bv .-1), (Nn.*Nm) ) .+ 1
    mv =  div.( (bv.-1 .-(lv.-1).*Nn.*Nm), Nn ) .+ 1
    nv =  bv.-(mv.-1).*Nn-(lv.-1).*Nn.*Nm
    ##  Gs have different shape than Us !!

    ####===================================

    ####================================================
    ##     Serial version
    ####================================================
    postm  = Array{Float64}(undef,Na)
    ddiff  = Array{Float64}(undef,Nb)
    Zh     = Array{Float64}(undef,Na)
    elUDZh = Array{Float64}(undef,Na)

    everynit = Na>20 ? div(Na,20) : 1

    startt = time()
    ## dobs - dcalc(mprior)
    @inbounds for b=1:Nb
        if (b%everynit==0) | (b==2)
            eta =  ( (time()-startt)/float(b-1) * (Na-b+1) ) /60.0
            reta = round(eta,digits=3)
            print("posteriormean(): loop 1/3, $b of $Nb; ETA: $reta min  \r")
            flush(stdout)
        end

        # !!---------------------------------------------------------------------------
        # ddiff(b) = dobs(b) - sum(mprior * G1(lv(b),iv) * G2(mv(b),jv) * G3(nv(b),kv))
        # !!---------------------------------------------------------------------------

        datp = 0.0
        @inbounds for j=1:Na
            elG = G1[lv[b],iv[j]] * G2[mv[b],jv[j]] * G3[nv[b],kv[j]]
            datp = datp +  mprior[j] * elG
        end
        ddiff[b] = dobs[b] - datp
    end

    startt = time()
    #tmpzhi = Array{Float64,1}(undef,Nb)
    @inbounds for i=1:Na

        if (i%everynit==0) | (i==2)
            eta =  ( (time()-startt)/float(i-1) * (Na-i+1) ) /60.0
            reta = round(eta,digits=3)
            print("posteriormean(): loop 2/3, $i of $Na; ETA: $reta min  \r")
            flush(stdout)
        end

        ## compute Zh
        # tmpzhi .= ddiff .* (Z1[iv[i],lv] .* Z2[jv[i],mv] .* Z3[kv[i],nv])
        # Zh[i] =  sum(tmpzhi)

        Zh[i]=0.0
        @inbounds  for j=1:Nb
            tZZ = Z1[iv[i],lv[j]] * Z2[jv[i],mv[j]] * Z3[kv[i],nv[j]]
            Zh[i] = Zh[i] + tZZ * ddiff[j]
        end
    end

    startt = time()
    ### need to re-loop because full Zh is needed
    @inbounds for i=1:Na

        if (i%everynit==0) | (i==2)
            eta =  ( (time()-startt)/float(i+1) * (Na-i+1) ) /60.0
            reta = round(eta,digits=3)
            print("posteriormean(): loop 3/3, $i of $Na; ETA: $reta min   \r")
            flush(stdout)
        end

        ## UD times Zh
        elUDZh[i] = 0.0
        @inbounds for j=1:Na
            # element of row of UD
            elrowUD = U1[iv[i],iv[j]] * U2[jv[i],jv[j]] *
                U3[kv[i],kv[j]] * diaginvlambda[j]
            # element of final vector
            elUDZh[i] = elUDZh[i] + elrowUD * Zh[j]
        end

        ## element of the posterior mean
        postm[i] = mprior[i] + elUDZh[i] # sum(bigmatrow.*ddiff)
    end
    println()

    return postm
end

##============================================================================

@doc raw"""
    blockpostcov_serial(klifac::KLIFactors,astart::Int64,aend::Int64,
                 bstart::Int64,bend::Int64 )

Computes a block of the posterior covariance, purely serial version.

# Arguments
- `klifac`: a structure containing the required "factors" previously computed with
    the function `calcfactors()`. It includes
    * U1,U2,U3 `` \mathbf{U}_1``, ``\mathbf{U}_2``, ``\mathbf{U}_3`` of ``F_{\sf{A}}``
    * diaginvlambda ``F_{\sf{B}}``
    * iUCm1, iUCm2, iUCm3  ``\mathbf{U}_1^{-1} \mathbf{C}_{\rm{M}}^{\rm{x}} ``,
      ``\mathbf{U}_2^{-1}  \mathbf{C}_{\rm{M}}^{\rm{y}}``,
      ``\mathbf{U}_2^{-1}  \mathbf{C}_{\rm{M}}^{\rm{z}}`` of  ``F_{\sf{C}} ``
- `Gfwd`: a structure containing the three forward model matrices  G1,G2,G3, where
     `` \mathbf{G} =  \mathbf{G_1} \otimes \mathbf{G_2} \otimes \mathbf{G_3} ``
- `astart, aend`: indices of the first and last rowa of the requested block
- `bstart, bend`: indices of the first and last columns of the requested block

# Returns
- The requested block of the posterior covariance.

"""
function blockpostcov_serial(klifac::KLIFactors,
                      astart::Int64,aend::Int64,
                      bstart::Int64,bend::Int64 )

    println("blockpostcov_serial(): Serial version")

    ##--------------
    U1,U2,U3 = klifac.U1,klifac.U2,klifac.U3
    diaginvlambda = klifac.invlambda
    iUCm1,iUCm2,iUCm3 = klifac.iUCm1,klifac.iUCm2,klifac.iUCm3

    ##-----------------
    Ni = size(U1,1)
    Nl = size(U1,2)
    Nj = size(U2,1)
    Nm = size(U2,2)
    Nk = size(U3,1)
    Nn = size(U3,2)
    Na = Ni*Nj*Nk
    Nb = Nl*Nm*Nn

    ## check limits of requested block
    if astart<1 || aend>Na || astart>aend || bstart<1 || bend>Na || bstart>bend
        error("blockpostcov(): Wrong size of the requested block array.")
    end

    ## sizes----
    # Nr12 = size(U1,2)*size(U2,2)*size(U3,2)
    # Nc1  = size(Z1,1)*size(Z2,1)*size(Z3,1)

    ##-----------------
    av = collect(1:Na)
    ## vectors containing all possible indices for
    ##    row calculations of Kron prod AxBxC
    iv =  div.( (av.-1), (Nk.*Nj) ) .+1
    jv =  div.( (av.-1 .-(iv.-1).*Nk.*Nj), Nk ) .+ 1
    kv =  av.-(jv.-1).*Nk.-(iv.-1).*Nk.*Nj

    nci = aend-astart+1
    ncj = bend-bstart+1
    postC = Array{Float64}(undef,nci,ncj)

    ####================================================
    ##     Serial version
    ####================================================
    row2  = Array{Float64}(undef,Na)

    @inbounds for a=astart:aend
        if a%100==0
            print("blockpostcov():  $a of $(astart) to $(aend) \r")
            flush(stdout)
        end
        ## calculate one row of first two factors
        ## row of  Kron prod AxBxC times a diag matrix (fb)
        @inbounds for q=1:Na
            row2[q] = U1[iv[a],iv[q]] * U2[jv[a],jv[q]] * U3[kv[a],kv[q]] * diaginvlambda[q]
        end
        @inbounds for b=bstart:bend
            postC[a,b] = 0.0
            @inbounds for p=1:Na
                ## calculate one column of fc
                col1 = iUCm1[iv[p],iv[b]] * iUCm2[jv[p],jv[b]] * iUCm3[kv[p],kv[b]]
                ## calculate one element
                postC[a,b] = postC[a,b] + row2[p] * col1
            end
        end
    end

    println()
    return postC
end

##==========================================================


text:
const WEATHER_DECADES = Dict(Date(1960) => "1960-1969",
                             Date(1970) => "1970-1979",
                             Date(1980) => "1980-1989",
                             Date(1990) => "1990-1999",
                             Date(2000) => "2000-2009",
                             Date(2010) => "2010-2018")

# Interface methods

rasterlayers(::Type{WorldClim{Weather}}) = (:tmin, :tmax, :prec)

function download_raster(T::Type{WorldClim{Weather}}, layers::Tuple=rasterlayers(T); kwargs...)
    all(l -> l in rasterlayers(T), layers) || throw(ArgumentError("Layers must be from $(rasterlayers(T))"))
    map(l -> download_raster(T, l; kwargs...), layers)
end
function download_raster(T::Type{WorldClim{Weather}}, layer::Symbol; dates)
    dates = _date_sequence(dates, Month(1))
    decadestart = Date.(1960:10:2020)
    raster_paths = String[]

    for i in eachindex(decadestart[1:end-1])
        # At least one date is in the decade
        any(d -> d >= decadestart[i] && d < decadestart[i+1], dates) || continue
        zip_path = zippath(T, layer, decadestart[i])
        _maybe_download(zipurl(T, layer, decadestart[i]), zip_path)
        zf = ZipFile.Reader(zip_path)
        for d in dates
            raster_path = rasterpath(T, layer, d)
            mkpath(dirname(raster_path))
            if !isfile(raster_path)
                raster_name = rastername(T, layer, d)
                println("Writing $(raster_path)...")
                write(raster_path, read(_zipfile_to_read(raster_name, zf)))
            end
            push!(raster_paths, raster_path)
        end
        close(zf)
    end
    return raster_paths
end

rastername(T::Type{<:WorldClim{Weather}}, key, date) =
    joinpath("wc2.1_2.5m_$(key)_$(_date2string(T, date)).tif")

zipname(T::Type{<:WorldClim{Weather}}, key, decade) =
    "wc2.1_2.5m_$(key)_$(WEATHER_DECADES[decade]).zip"

zipurl(T::Type{<:WorldClim{Weather}}, key, decade) =
    joinpath(WORLDCLIM_URL, "hist", zipname(T, key, decade))

zippath(T::Type{<:WorldClim{Weather}}, key, decade) =
    joinpath(rasterpath(T), "zips", zipname(T, key, decade))

# Utility methods

_dateformat(::Type{<:WorldClim}) = DateFormat("yyyy-mm")

_filename2date(T::Type{<:WorldClim}, fn::AbstractString) =
    _string2date(T, basename(fn)[findfirst(r"\d\d\d\d-\d\d", basename(fn))])


text:
@testset "KLIEP -- $optlib" for optlib in [OptimLib, ConvexLib]
  for (i, (pair, rtol)) in enumerate([(pair₁, 2e-1), (pair₂, 4e-1)])
    d_nu, d_de = pair
    rng = MersenneTwister(42)
    x_nu, x_de = rand(rng, d_nu, 1_000), rand(rng, d_de, 500)

    # estimated density ratio
    kliep = KLIEP(σ=1.0, b=100, rng=rng)
    r̂ = densratio(x_nu, x_de, kliep, optlib=optlib)

    # simplex constraints
    @test abs(mean(r̂) - 1) ≤ 1e-2
    @test all(r̂ .≤ Inf)

    r = pdf.(d_nu, x_de) ./ pdf.(d_de, x_de)
    @test r ≈ r̂ rtol=rtol

    if visualtests
      gr(size=(800, 800))
      @test_reference "data/KLIEP-$optlib-$i.png" plot_d_nu(pair, x_de, r̂)
    end
  end
end


text:
###=============================================================================
#
#          FILE: CloudArray.jl
#
#         USAGE: include("CloudArray.jl")
#
#   DESCRIPTION: Implements the CloudArrays functionalities for Julia
#
#       OPTIONS: ---
#  DEPENDENCIES: DistributedArrays
#          BUGS: ---
#         NOTES: ---
#        AUTHOR: Naelson Douglas C. Oliveira
#  ORGANIZATION: GSD-UFAL
#       CREATED: 2015-10-04 16:32
###=============================================================================

###=========================   Environment setup   =============================

#set this as true if you want to get the execution time of the CArray (the default is false), if true our constructors outputs will be
#data, number of chunks created, time to create the container   eg.: darray,n_containers, time = DArray("floats.txt")
#if false it will return only the DArray

module CloudArray

    using DistributedArrays

    include("Infra.jl")

    #You can set this var as true to test the cloudarray without the Infra.jl management. Keep in mind that instead of using the cloud to share your data, it will use localprocs with addprocs(1). This is only for tests purposes.
    local_workers = false

    ###=============================================================================
    ###=========================   Auxiliary functions   ===========================

    @doc """
    ### rmall()
    Remove all workers.
    The same as rmprocs(workers())
    WARNING: Every data stored in the workers will be lost.
    """ ->
    #Removes all the workers()
    function rmprocs_all()
        rmprocs(workers())
    end

    @doc """
    ### compose(carray1::DArray, carray2::DArray)
    Gets 2 DArrays and return both as one. No data IO is involved, only the references
    """ ->
    function compose(array1::DistributedArrays.DArray, carray2::DistributedArrays.DArray)
        return DArray([carray1.chunks;carray2.chunks])
    end

    @doc """
    ### lines_producer(input::String)
    Takes an input filepath and returns a task which generetas one of the file by call.
    ## use consume(::Task) to get the lines in the task
    """->
    function lines_producer(input::AbstractString="floats.txt")
        file = open(input)
        for line in eachline(file)
            produce(line)
        end
        close(file)
    end

    @doc """
    Gets two values, low and high and produces these values one by one. If low > high, the are swaped
    """->
    function seq_producer(low::Int, high::Int)

        if (low > high)
            low, high = high, low
        end

        for i=low:high
            produce(i)
        end
    end

    @doc """
    An auxiliary function for the CloudArray. It gets a line of data, which was previously generated by a task, and parse this line as a float or only removes the eol's from the line and outputs it as a string.
    """->
    function parse_by_type(data, is_numeric::Bool=false)
        if (is_numeric)
            return 	parse(chomp(string(data)))
        else
            return chomp(string(data))
        end
    end

    @doc """
    Gets an Array produces it with produce(value) for each line of the Array.
    This function is use to cread tasks in the task_from_array(...) function
    """->
    function array_producer(a::AbstractArray)
               for i in a
                   produce(i)
           end
    end

    @doc """
    ### task_from_text(input::AbstractString)
    Gets a filepath as input and returns a task which produces the lines from the file.
    """->
    function task_from_text(input::AbstractString="floats.txt",delim::Char='\n')
        return ( Task( () -> lines_producer(input)))
    end

    @doc """
    Gets  an Array as input and creates a task whose will produce every values in the array one by one with consume(::Task) calls
    """->
    function task_from_array(input::AbstractArray)
        return(Task(()->array_producer(input)))
    end

    @doc """
    The same as DistributedArrays.distribute(::AbstractArray), but using the CloudArray management.
    """->
    distribute_cloud(input::AbstractArray) = carray_from_task(task_from_array(input))

    @doc """
    Creates a CloudArray from an Array. It's actualt a call of distribute_cloud(input::AbstractArray) with other name
    """->
    DistributedArrays.DArray(input::AbstractArray) = begin distribute_cloud(input) end

    ###=============================================================================
    ###=========================   Core Constructors   =============================

    @doc """
    Core constructor of the CloudArray

    Creates a CloudArray from a task. The task must generate data (lines) per consume(::Task) call. Each call will be stored as a CloudArray element.

    * `is_numeric::Bool` should be set to `true` if your data is composed by numbers. The default value is `true`.
    * `chunk_max_size::Int` means the maximum size of each CloudArray chunk in bytes. The default value is 1MB (1024*1024).
    * `debug::Bool` enables debug mode, default is `false`. If you set `true`, the function output will be `DArray, [number_of_chunks], [time_to_instance_each_chunk]` instead of only returning a `DArray`.

    This constructor was created with tasks to make easier to tune a custom CloudArray. The user just needs to create a task which generates the lines of the CloudArray by each call, e.g.: `task_from_array(...)` and `task_from_text(...)`.
    """->
    function carray_from_task(generator::Task=task_from_text("floats.txt"), chunk_max_size::Int=1024*1024; is_numeric::Bool=true, debug::Bool=false)

        plot_data =  []
        containers = []
        time = []
        refs = Array(RemoteRef,0)
        while (generator.state != :done) #while there's something to be read in the source

            #creates an array to buffer the lines from the source
            if (is_numeric)
                local_array = Array(Float64,0)
            else
                local_array = Array(AbstractString,0)
            end

            while ((sizeof(local_array) < chunk_max_size) & (generator.state != :done)) #while the current container still have free space, store the lines in the buffer

                line = consume(generator)
                if (typeof(line) != Void && line !='\n') #we don't want empty lines in our CArray, right?

                push!(local_array,parse_by_type(line,is_numeric))
                end

            end
            size_local = length(local_array)

            tic()
            if !local_workers

                create_containers(1,0,250,tunnel=true)
                #push!(plot_data, [time;workers])
            else
                addprocs(1)
            end
            println("New worker added\nTotal: ",nworkers())

            push!(containers,ncontainers())
            push!(time,toc())
            @everywhere using DistributedArrays

            container_chunk = distribute(local_array, procs=[workers()[nworkers()]])  #stores a chunk(as an independent darray) in a container
            refs = vcat(refs, container_chunk.chunks)

        end

        if debug
            return DArray(refs),containers,time #this return is used for debugs and tests purposes.
        end
            return DArray(refs) #This is the standard return
    end

    ###=============================================================================
    ###=========================   Constructors   ==================================

    @doc """
    ### DArray(input::AbstractString,args...)
    creates a CloudArray from the file used on the input
    """ ->
    DistributedArrays.DArray(input::AbstractString, args...) = begin
        return	carray_from_task(task_from_text(input), args...)
    end

    DistributedArrays.DArray(input::Array, args...) = begin
        return	carray_from_task(task_from_array(input), is_numeric=true, args...)
    end

    ###=============================================================================

end


text:
numeric_type(::Interval{T}) where {T} = T
numeric_type(::IntervalBox{N, T}) where {N, T} = T

"""
    minimise(f, X, structure = SortedVector, tol=1e-3)
    or
    minimise(f, X, structure = HeapedVector, tol=1e-3)
    or
    minimise(f, X, tol=1e-3) in this case the default value of "structure" is "HeapedVector"

Find the global minimum of the function `f` over the `Interval` or `IntervalBox` `X`
using the Moore-Skelboe algorithm. The way in which vector elements are
kept arranged can be either a heaped array or a sorted array.
If you not specify any particular strategy to keep vector elements arranged then
by default heaped array is used.

For higher-dimensional functions ``f:\\mathbb{R}^n \\to \\mathbb{R}``, `f` must take a single vector argument.

Returns an interval containing the global minimum, and a list of boxes that contain the minimisers.
"""
function minimise(f, X::T; structure = HeapedVector, tol=1e-3) where {T}
    nT = numeric_type(X)

    # list of boxes with corresponding lower bound, arranged according to selected structure :
    working = structure([(X, nT(∞))], x->x[2])
    minimizers = T[]
    global_min = nT(∞)  # upper bound

    num_bisections = 0

    while !isempty(working)

        (X, X_min) = popfirst!(working)

        if X_min > global_min    # X_min == inf(f(X))
            continue
        end

        # find candidate for upper bound of global minimum by just evaluating a point in the interval:
        m = sup(f(Interval.(mid.(X))))   # evaluate at midpoint of current interval

        if m < global_min
            global_min = m
        end

        # Remove all boxes whose lower bound is greater than the current one:
        filter_elements!(working , (X, global_min) )

        if diam(X) < tol
            push!(minimizers, X)
        else
            X1, X2 = bisect(X)
            push!( working, (X1, inf(f(X1))) )
            push!( working, (X2, inf(f(X2))) )
            num_bisections += 1
        end

    end

    lower_bound = minimum(inf.(f.(minimizers)))

    return Interval(lower_bound, global_min), minimizers
end

"""
    maximise(f, X, structure = SortedVector, tol=1e-3)
    or
    maximise(f, X, structure = HeapedVector, tol=1e-3)
    or
    maximise(f, X, tol=1e-3) in this case the default value of "structure" is "HeapedVector"

Find the global maximum of the function `f` over the `Interval` or `IntervalBox` `X`
using the Moore-Skelboe algorithm. See [`minimise`](@ref) for a description
of the available options.
"""
function maximise(f, X::T; structure=HeapedVector, tol=1e-3) where {T}
    bound, minimizers = minimise(x -> -f(x), X, structure=structure, tol=tol)
    return -bound, minimizers
end


text:
# ------------------------------------------------------------------
# Licensed under the ISC License. See LICENSE in the project root.
# ------------------------------------------------------------------

"""
    AbstractRegion

A region in a `N`-dimensional space with coordinates of type `T`.
"""
abstract type AbstractRegion{T,N} end

#------------------
# IMPLEMENTATIONS
#------------------
include("regions/rectangle.jl")


text:
# Object definition
struct Coupler{N, N1}
    a::NTuple{N,Float64}
    b::NTuple{N1,Float64}
    Z₀::Float64
    function Coupler(a, b, Z₀)
        N = length(a)
        N1 = length(b)
        N == N1 + 1 ||
            error("Vector of series transmission lines must be 1 shorter than vector of branches")
        new{N,N1}(a, b, Z₀)
    end
end

# Omitting the impedance
Coupler(a::NTuple, b::NTuple) = Coupler(a, b, 1.0)

# Omitting b
Coupler(a, Z₀::T where {T<:Real} = 1.0) = Coupler(a, ntuple(Returns(Z₀),length(a)-1), Z₀)

Coupler(a::AbstractArray{<:Real}, b::T where{T<:Tuple}, x...) = Coupler(tuple(a...), b, x...)
Coupler(a, b::AbstractArray{<:Real}, x...) = Coupler(a, tuple(b...), x...)

a = [0.0965767, 0.114064, 0.131845, 0.114064, 0.0965767]
b  = [1.02556, 1.04526, 1.04526, 1.02556]
c = Coupler(a,b)


text:
using Logging

const date_re = r"(\w{3}\s\w{3}\s{1,2}\d{1,2}\s\d{2}:\d{2}:\d{2}\s\d{4})"

const file_re = r"\s{2}"*         # file lines start with 2 spaces
    r"(.*?)\s+"*            # capture filename non-greedy, eating remaining spaces
    r"([ADHSR]*)"*          # capture file mode
    r"\s+"*                 # after the mode you can have any number of spaces
    r"(\d+)"*               # file size
    r"\s+"*                 # spaces after file size
    date_re*
    r"$"

const path_dir_splitter = r"^(.*?)([/\\]+)([^/\\]*)$"

mutable struct SmbPath
    host::String
    auth_file::String
    path::String
end

function run_smbclient(path::SmbPath, command::String; dir=nothing)
    if dir === nothing
        cmd = `smbclient -A $(path.auth_file) $(path.host) -c $command`
    else
        cmd = `smbclient -A $(path.auth_file) $(path.host) -D $dir  -c $command`
    end
    @info "smbclient" cmd
    output = readlines(cmd)
end

function Base.Filesystem.readdir(path::SmbPath)
    output = run_smbclient(path, "ls", dir=path.path)
    result = String[]
    for l in output
        res = match(file_re, l)
        if res !== nothing && res[1] != "." && res[1] != ".."
            push!(result, res[1])
        else
            nothing
        end
    end
    result
end

function allinfo(path::SmbPath)
    smb_path = path.path
    output = run_smbclient(path, "allinfo $smb_path")
    infodict = Dict{String, Any}()
    for l in output
        attr_re = r"attributes:\s+([ADHSR]*)\s+\([0-9]+\)"
        res = match(attr_re, l)
        if res !== nothing
            infodict["attributes"] = res[1]
        end
    end
    infodict
end

function Base.Filesystem.isdir(path::SmbPath)
    info = allinfo(path)
    if haskey(info, "attributes")
        'D' in info["attributes"]
    else
        false
    end
end

function Base.Filesystem.mkdir(path::SmbPath, name::String)
    output = run_smbclient(path, "mkdir $name", dir=path.path)
end

function Base.Filesystem.rm(path::SmbPath; recursive::Bool=false)
    if recursive
        for (root, dirs, files) in walkdir(path, topdown=false)
            for dir in dirs
                rm(joinpath(root, dir))
            end
            for file in files
                rm(joinpath(root, file))
            end
        end
    end

    dir, name = splitdir(path)
    if isdir(path)
        command = "rmdir $name"
        if !isempty(readdir(path))
            error("directory is not empty")
        end
    else
        command = "rm $name"
    end
    output = run_smbclient(path, command, dir=dir.path)
end

function Base.Filesystem.mkdir(path::SmbPath)
    dname, basename = Base.Filesystem.splitdir(path)
    mkdir(dname, basename)
end

function Base.Filesystem.cp(src::String, dst::SmbPath)
    dir, name = splitdir(dst)
    run_smbclient(dst, "put $src $name", dir=dir.path)
end

function Base.Filesystem.joinpath(root::SmbPath, parts::AbstractString...)
    comps = [root.path, parts...]
    SmbPath(root.host, root.auth_file, join(comps, '\\'))
end

function splitdir_str(a::String, b::String)
    m = match(path_dir_splitter,b)
    m === nothing && return (a,b)
    a = string(a, isempty(m.captures[1]) ? m.captures[2][1] : m.captures[1])
    a, String(m.captures[3])
end

function Base.Filesystem.splitdir(path::SmbPath)
    (a, b) = splitdir_str("", path.path)
    SmbPath(path.host, path.auth_file, a), b
end

function Base.Filesystem.dirname(path::SmbPath)
    dirn = splitdir_str("", path.path)[1]
    SmbPath(path.host, path.auth_file, dirn)
end

function Base.Filesystem.basename(path::SmbPath)
    splitdir(path.path)[2]
end

function Base.Filesystem.mkpath(path::SmbPath)
    if !isdir(path)
        dname = Base.Filesystem.dirname(path)
        mkpath(dname)
        mkdir(path)
    end
end

Base.Filesystem.islink(path::SmbPath) = false

Base.Filesystem.isabspath(path::SmbPath) = true

host = "//sun.amolf.nl/tenwolde"
sun_home = SmbPath(host, "/home/ipausers/reinhardt/sun_auth.txt", "home-folder\\reinhardt")


text:
# This file is auto-generated by AWSMetadata.jl
using AWS
using AWS.AWSServices: kinesis_video
using AWS.Compat
using AWS.UUIDs

"""
    CreateSignalingChannel()

Creates a signaling channel.   CreateSignalingChannel is an asynchronous operation.

# Required Parameters
- `ChannelName`: A name for the signaling channel that you are creating. It must be unique for each AWS account and AWS Region.

# Optional Parameters
- `ChannelType`: A type of the signaling channel that you are creating. Currently, SINGLE_MASTER is the only supported channel type.
- `SingleMasterConfiguration`: A structure containing the configuration for the SINGLE_MASTER channel type.
- `Tags`: A set of tags (key-value pairs) that you want to associate with this channel.
"""
create_signaling_channel(ChannelName; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/createSignalingChannel", Dict{String, Any}("ChannelName"=>ChannelName); aws_config=aws_config)
create_signaling_channel(ChannelName, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/createSignalingChannel", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ChannelName"=>ChannelName), args)); aws_config=aws_config)

"""
    CreateStream()

Creates a new Kinesis video stream.  When you create a new stream, Kinesis Video Streams assigns it a version number. When you change the stream's metadata, Kinesis Video Streams updates the version.   CreateStream is an asynchronous operation. For information about how the service works, see How it Works.  You must have permissions for the KinesisVideo:CreateStream action.

# Required Parameters
- `StreamName`: A name for the stream that you are creating. The stream name is an identifier for the stream, and must be unique for each account and region.

# Optional Parameters
- `DataRetentionInHours`: The number of hours that you want to retain the data in the stream. Kinesis Video Streams retains the data in a data store that is associated with the stream. The default value is 0, indicating that the stream does not persist data. When the DataRetentionInHours value is 0, consumers can still consume the fragments that remain in the service host buffer, which has a retention time limit of 5 minutes and a retention memory limit of 200 MB. Fragments are removed from the buffer when either limit is reached.
- `DeviceName`: The name of the device that is writing to the stream.   In the current implementation, Kinesis Video Streams does not use this name.
- `KmsKeyId`: The ID of the AWS Key Management Service (AWS KMS) key that you want Kinesis Video Streams to use to encrypt stream data. If no key ID is specified, the default, Kinesis Video-managed key (aws/kinesisvideo) is used.  For more information, see DescribeKey.
- `MediaType`: The media type of the stream. Consumers of the stream can use this information when processing the stream. For more information about media types, see Media Types. If you choose to specify the MediaType, see Naming Requirements for guidelines. Example valid values include \"video/h264\" and \"video/h264,audio/aac\". This parameter is optional; the default value is null (or empty in JSON).
- `Tags`: A list of tags to associate with the specified stream. Each tag is a key-value pair (the value is optional).
"""
create_stream(StreamName; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/createStream", Dict{String, Any}("StreamName"=>StreamName); aws_config=aws_config)
create_stream(StreamName, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/createStream", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("StreamName"=>StreamName), args)); aws_config=aws_config)

"""
    DeleteSignalingChannel()

Deletes a specified signaling channel. DeleteSignalingChannel is an asynchronous operation. If you don't specify the channel's current version, the most recent version is deleted.

# Required Parameters
- `ChannelARN`: The Amazon Resource Name (ARN) of the signaling channel that you want to delete.

# Optional Parameters
- `CurrentVersion`: The current version of the signaling channel that you want to delete. You can obtain the current version by invoking the DescribeSignalingChannel or ListSignalingChannels API operations.
"""
delete_signaling_channel(ChannelARN; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/deleteSignalingChannel", Dict{String, Any}("ChannelARN"=>ChannelARN); aws_config=aws_config)
delete_signaling_channel(ChannelARN, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/deleteSignalingChannel", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ChannelARN"=>ChannelARN), args)); aws_config=aws_config)

"""
    DeleteStream()

Deletes a Kinesis video stream and the data contained in the stream.  This method marks the stream for deletion, and makes the data in the stream inaccessible immediately.    To ensure that you have the latest version of the stream before deleting it, you can specify the stream version. Kinesis Video Streams assigns a version to each stream. When you update a stream, Kinesis Video Streams assigns a new version number. To get the latest stream version, use the DescribeStream API.  This operation requires permission for the KinesisVideo:DeleteStream action.

# Required Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the stream that you want to delete.

# Optional Parameters
- `CurrentVersion`: Optional: The version of the stream that you want to delete.  Specify the version as a safeguard to ensure that your are deleting the correct stream. To get the stream version, use the DescribeStream API. If not specified, only the CreationTime is checked before deleting the stream.
"""
delete_stream(StreamARN; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/deleteStream", Dict{String, Any}("StreamARN"=>StreamARN); aws_config=aws_config)
delete_stream(StreamARN, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/deleteStream", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("StreamARN"=>StreamARN), args)); aws_config=aws_config)

"""
    DescribeSignalingChannel()

Returns the most current information about the signaling channel. You must specify either the name or the Amazon Resource Name (ARN) of the channel that you want to describe.

# Optional Parameters
- `ChannelARN`: The ARN of the signaling channel that you want to describe.
- `ChannelName`: The name of the signaling channel that you want to describe.
"""
describe_signaling_channel(; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/describeSignalingChannel"; aws_config=aws_config)
describe_signaling_channel(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/describeSignalingChannel", args; aws_config=aws_config)

"""
    DescribeStream()

Returns the most current information about the specified stream. You must specify either the StreamName or the StreamARN.

# Optional Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the stream.
- `StreamName`: The name of the stream.
"""
describe_stream(; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/describeStream"; aws_config=aws_config)
describe_stream(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/describeStream", args; aws_config=aws_config)

"""
    GetDataEndpoint()

Gets an endpoint for a specified stream for either reading or writing. Use this endpoint in your application to read from the specified stream (using the GetMedia or GetMediaForFragmentList operations) or write to it (using the PutMedia operation).   The returned endpoint does not have the API name appended. The client needs to add the API name to the returned endpoint.  In the request, specify the stream either by StreamName or StreamARN.

# Required Parameters
- `APIName`: The name of the API action for which to get an endpoint.

# Optional Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the stream that you want to get the endpoint for. You must specify either this parameter or a StreamName in the request.
- `StreamName`: The name of the stream that you want to get the endpoint for. You must specify either this parameter or a StreamARN in the request.
"""
get_data_endpoint(APIName; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/getDataEndpoint", Dict{String, Any}("APIName"=>APIName); aws_config=aws_config)
get_data_endpoint(APIName, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/getDataEndpoint", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("APIName"=>APIName), args)); aws_config=aws_config)

"""
    GetSignalingChannelEndpoint()

Provides an endpoint for the specified signaling channel to send and receive messages. This API uses the SingleMasterChannelEndpointConfiguration input parameter, which consists of the Protocols and Role properties.  Protocols is used to determine the communication mechanism. For example, if you specify WSS as the protocol, this API produces a secure websocket endpoint. If you specify HTTPS as the protocol, this API generates an HTTPS endpoint.   Role determines the messaging permissions. A MASTER role results in this API generating an endpoint that a client can use to communicate with any of the viewers on the channel. A VIEWER role results in this API generating an endpoint that a client can use to communicate only with a MASTER.

# Required Parameters
- `ChannelARN`: The Amazon Resource Name (ARN) of the signalling channel for which you want to get an endpoint.

# Optional Parameters
- `SingleMasterChannelEndpointConfiguration`: A structure containing the endpoint configuration for the SINGLE_MASTER channel type.
"""
get_signaling_channel_endpoint(ChannelARN; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/getSignalingChannelEndpoint", Dict{String, Any}("ChannelARN"=>ChannelARN); aws_config=aws_config)
get_signaling_channel_endpoint(ChannelARN, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/getSignalingChannelEndpoint", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ChannelARN"=>ChannelARN), args)); aws_config=aws_config)

"""
    ListSignalingChannels()

Returns an array of ChannelInfo objects. Each object describes a signaling channel. To retrieve only those channels that satisfy a specific condition, you can specify a ChannelNameCondition.

# Optional Parameters
- `ChannelNameCondition`: Optional: Returns only the channels that satisfy a specific condition.
- `MaxResults`: The maximum number of channels to return in the response. The default is 500.
- `NextToken`: If you specify this parameter, when the result of a ListSignalingChannels operation is truncated, the call returns the NextToken in the response. To get another batch of channels, provide this token in your next request.
"""
list_signaling_channels(; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listSignalingChannels"; aws_config=aws_config)
list_signaling_channels(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listSignalingChannels", args; aws_config=aws_config)

"""
    ListStreams()

Returns an array of StreamInfo objects. Each object describes a stream. To retrieve only streams that satisfy a specific condition, you can specify a StreamNameCondition.

# Optional Parameters
- `MaxResults`: The maximum number of streams to return in the response. The default is 10,000.
- `NextToken`: If you specify this parameter, when the result of a ListStreams operation is truncated, the call returns the NextToken in the response. To get another batch of streams, provide this token in your next request.
- `StreamNameCondition`: Optional: Returns only streams that satisfy a specific condition. Currently, you can specify only the prefix of a stream name as a condition.
"""
list_streams(; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listStreams"; aws_config=aws_config)
list_streams(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listStreams", args; aws_config=aws_config)

"""
    ListTagsForResource()

Returns a list of tags associated with the specified signaling channel.

# Required Parameters
- `ResourceARN`: The Amazon Resource Name (ARN) of the signaling channel for which you want to list tags.

# Optional Parameters
- `NextToken`: If you specify this parameter and the result of a ListTagsForResource call is truncated, the response includes a token that you can use in the next request to fetch the next batch of tags.
"""
list_tags_for_resource(ResourceARN; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/ListTagsForResource", Dict{String, Any}("ResourceARN"=>ResourceARN); aws_config=aws_config)
list_tags_for_resource(ResourceARN, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/ListTagsForResource", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ResourceARN"=>ResourceARN), args)); aws_config=aws_config)

"""
    ListTagsForStream()

Returns a list of tags associated with the specified stream. In the request, you must specify either the StreamName or the StreamARN.

# Optional Parameters
- `NextToken`: If you specify this parameter and the result of a ListTagsForStream call is truncated, the response includes a token that you can use in the next request to fetch the next batch of tags.
- `StreamARN`: The Amazon Resource Name (ARN) of the stream that you want to list tags for.
- `StreamName`: The name of the stream that you want to list tags for.
"""
list_tags_for_stream(; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listTagsForStream"; aws_config=aws_config)
list_tags_for_stream(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/listTagsForStream", args; aws_config=aws_config)

"""
    TagResource()

Adds one or more tags to a signaling channel. A tag is a key-value pair (the value is optional) that you can define and assign to AWS resources. If you specify a tag that already exists, the tag value is replaced with the value that you specify in the request. For more information, see Using Cost Allocation Tags in the AWS Billing and Cost Management User Guide.

# Required Parameters
- `ResourceARN`: The Amazon Resource Name (ARN) of the signaling channel to which you want to add tags.
- `Tags`: A list of tags to associate with the specified signaling channel. Each tag is a key-value pair.

"""
tag_resource(ResourceARN, Tags; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/TagResource", Dict{String, Any}("ResourceARN"=>ResourceARN, "Tags"=>Tags); aws_config=aws_config)
tag_resource(ResourceARN, Tags, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/TagResource", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ResourceARN"=>ResourceARN, "Tags"=>Tags), args)); aws_config=aws_config)

"""
    TagStream()

Adds one or more tags to a stream. A tag is a key-value pair (the value is optional) that you can define and assign to AWS resources. If you specify a tag that already exists, the tag value is replaced with the value that you specify in the request. For more information, see Using Cost Allocation Tags in the AWS Billing and Cost Management User Guide.  You must provide either the StreamName or the StreamARN. This operation requires permission for the KinesisVideo:TagStream action. Kinesis video streams support up to 50 tags.

# Required Parameters
- `Tags`: A list of tags to associate with the specified stream. Each tag is a key-value pair (the value is optional).

# Optional Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the resource that you want to add the tag or tags to.
- `StreamName`: The name of the stream that you want to add the tag or tags to.
"""
tag_stream(Tags; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/tagStream", Dict{String, Any}("Tags"=>Tags); aws_config=aws_config)
tag_stream(Tags, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/tagStream", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Tags"=>Tags), args)); aws_config=aws_config)

"""
    UntagResource()

Removes one or more tags from a signaling channel. In the request, specify only a tag key or keys; don't specify the value. If you specify a tag key that does not exist, it's ignored.

# Required Parameters
- `ResourceARN`: The Amazon Resource Name (ARN) of the signaling channel from which you want to remove tags.
- `TagKeyList`: A list of the keys of the tags that you want to remove.

"""
untag_resource(ResourceARN, TagKeyList; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/UntagResource", Dict{String, Any}("ResourceARN"=>ResourceARN, "TagKeyList"=>TagKeyList); aws_config=aws_config)
untag_resource(ResourceARN, TagKeyList, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/UntagResource", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ResourceARN"=>ResourceARN, "TagKeyList"=>TagKeyList), args)); aws_config=aws_config)

"""
    UntagStream()

Removes one or more tags from a stream. In the request, specify only a tag key or keys; don't specify the value. If you specify a tag key that does not exist, it's ignored. In the request, you must provide the StreamName or StreamARN.

# Required Parameters
- `TagKeyList`: A list of the keys of the tags that you want to remove.

# Optional Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the stream that you want to remove tags from.
- `StreamName`: The name of the stream that you want to remove tags from.
"""
untag_stream(TagKeyList; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/untagStream", Dict{String, Any}("TagKeyList"=>TagKeyList); aws_config=aws_config)
untag_stream(TagKeyList, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/untagStream", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("TagKeyList"=>TagKeyList), args)); aws_config=aws_config)

"""
    UpdateDataRetention()

 Increases or decreases the stream's data retention period by the value that you specify. To indicate whether you want to increase or decrease the data retention period, specify the Operation parameter in the request body. In the request, you must specify either the StreamName or the StreamARN.   The retention period that you specify replaces the current value.  This operation requires permission for the KinesisVideo:UpdateDataRetention action. Changing the data retention period affects the data in the stream as follows:   If the data retention period is increased, existing data is retained for the new retention period. For example, if the data retention period is increased from one hour to seven hours, all existing data is retained for seven hours.   If the data retention period is decreased, existing data is retained for the new retention period. For example, if the data retention period is decreased from seven hours to one hour, all existing data is retained for one hour, and any data older than one hour is deleted immediately.

# Required Parameters
- `CurrentVersion`: The version of the stream whose retention period you want to change. To get the version, call either the DescribeStream or the ListStreams API.
- `DataRetentionChangeInHours`: The retention period, in hours. The value you specify replaces the current value. The maximum value for this parameter is 87600 (ten years).
- `Operation`: Indicates whether you want to increase or decrease the retention period.

# Optional Parameters
- `StreamARN`: The Amazon Resource Name (ARN) of the stream whose retention period you want to change.
- `StreamName`: The name of the stream whose retention period you want to change.
"""
update_data_retention(CurrentVersion, DataRetentionChangeInHours, Operation; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateDataRetention", Dict{String, Any}("CurrentVersion"=>CurrentVersion, "DataRetentionChangeInHours"=>DataRetentionChangeInHours, "Operation"=>Operation); aws_config=aws_config)
update_data_retention(CurrentVersion, DataRetentionChangeInHours, Operation, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateDataRetention", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("CurrentVersion"=>CurrentVersion, "DataRetentionChangeInHours"=>DataRetentionChangeInHours, "Operation"=>Operation), args)); aws_config=aws_config)

"""
    UpdateSignalingChannel()

Updates the existing signaling channel. This is an asynchronous operation and takes time to complete.  If the MessageTtlSeconds value is updated (either increased or reduced), it only applies to new messages sent via this channel after it's been updated. Existing messages are still expired as per the previous MessageTtlSeconds value.

# Required Parameters
- `ChannelARN`: The Amazon Resource Name (ARN) of the signaling channel that you want to update.
- `CurrentVersion`: The current version of the signaling channel that you want to update.

# Optional Parameters
- `SingleMasterConfiguration`: The structure containing the configuration for the SINGLE_MASTER type of the signaling channel that you want to update.
"""
update_signaling_channel(ChannelARN, CurrentVersion; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateSignalingChannel", Dict{String, Any}("ChannelARN"=>ChannelARN, "CurrentVersion"=>CurrentVersion); aws_config=aws_config)
update_signaling_channel(ChannelARN, CurrentVersion, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateSignalingChannel", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("ChannelARN"=>ChannelARN, "CurrentVersion"=>CurrentVersion), args)); aws_config=aws_config)

"""
    UpdateStream()

Updates stream metadata, such as the device name and media type. You must provide the stream name or the Amazon Resource Name (ARN) of the stream. To make sure that you have the latest version of the stream before updating it, you can specify the stream version. Kinesis Video Streams assigns a version to each stream. When you update a stream, Kinesis Video Streams assigns a new version number. To get the latest stream version, use the DescribeStream API.   UpdateStream is an asynchronous operation, and takes time to complete.

# Required Parameters
- `CurrentVersion`: The version of the stream whose metadata you want to update.

# Optional Parameters
- `DeviceName`: The name of the device that is writing to the stream.    In the current implementation, Kinesis Video Streams does not use this name.
- `MediaType`: The stream's media type. Use MediaType to specify the type of content that the stream contains to the consumers of the stream. For more information about media types, see Media Types. If you choose to specify the MediaType, see Naming Requirements. To play video on the console, you must specify the correct video type. For example, if the video in the stream is H.264, specify video/h264 as the MediaType.
- `StreamARN`: The ARN of the stream whose metadata you want to update.
- `StreamName`: The name of the stream whose metadata you want to update. The stream name is an identifier for the stream, and must be unique for each account and region.
"""
update_stream(CurrentVersion; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateStream", Dict{String, Any}("CurrentVersion"=>CurrentVersion); aws_config=aws_config)
update_stream(CurrentVersion, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = kinesis_video("POST", "/updateStream", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("CurrentVersion"=>CurrentVersion), args)); aws_config=aws_config)


text:
function imagesc(A::AbstractMatrix)
    Ny, Nx = size(A)

    if Nx > 30
        k = ceil(Nx / 30)
        A = A[:, 1:k:Nx]
    end

    if Ny > 30
        k = ceil(Ny / 30)
        A = A[1:k:Ny, :]
    end

    Ny, Nx = size(A)

    io = IOBuffer()

    t1 = quantile(A[:], 0.2)
    t2 = quantile(A[:], 0.4)
    t3 = quantile(A[:], 0.6)
    t4 = quantile(A[:], 0.8)

    print(io, "\n\n")

    for i = 1:Ny
        print(io, "\t")
        for j = 1:Nx
            if A[i, j] <= t1
                print(io, "  ")
            elseif A[i, j] <= t2
                print(io, ". ")
            elseif A[i, j] <= t3
                print(io, "+ ")
            elseif A[i, j] <= t4
                print(io, "# ")
            else
                print(io, "@#")
            end
        end
        print(io, "\n")
    end

    print(io, "\n")

    return ASCIIPlot(bytestring(io))
end
function imagesc(A::AbstractMatrix)
    Ny, Nx = size(A)

    if Nx > 30
        k = ceil(Nx / 30)
        A = A[:, 1:k:Nx]
    end

    if Ny > 30
        k = ceil(Ny / 30)
        A = A[1:k:Ny, :]
    end

    Ny, Nx = size(A)

    io = IOBuffer()

    t1 = quantile(A[:], 0.2)
    t2 = quantile(A[:], 0.4)
    t3 = quantile(A[:], 0.6)
    t4 = quantile(A[:], 0.8)

    print(io, "\n\n")

    for i = 1:Ny
        print(io, "\t")
        for j = 1:Nx
            if A[i, j] <= t1
                print(io, "  ")
            elseif A[i, j] <= t2
                print(io, ". ")
            elseif A[i, j] <= t3
                print(io, "+ ")
            elseif A[i, j] <= t4
                print(io, "# ")
            else
                print(io, "@#")
            end
        end
        print(io, "\n")
    end

    print(io, "\n")

    return ASCIIPlot(bytestring(io))
end


text:
using Test
using TightBindingLattice

@testset "sitepermutation" begin
    @test_throws ArgumentError SitePermutation([4,3])
    @test_throws ArgumentError SitePermutation([4,3,2,2])
    @test isidentity(SitePermutation(Int[]))
    @test isidentity(SitePermutation([1,2,3]))
    @test !isidentity(SitePermutation([3,2,1]))

    op1 = SitePermutation([3,1,2]) # 1=>3, 2=>1, 3=>2
    op1p = SitePermutation([3,1,2])
    op2 = SitePermutation([2,1,3]) # 1=>2, 2=>1, 3=>3

    @test_throws BoundsError op1(10)
    @test op1(1) == 3 && op1(2) == 1 && op1(3) == 2

    @test op2 * op1 != SitePermutation([1,2,3])
    @test op2 * op1 == SitePermutation([3,2,1])

    @test op1 !== op1p
    @test op1 == op1p
    @test hash(op1) == hash(op1p)

    @test op1 != op2
    @test hash(op1) != hash(op2)

    @test inv(op1) == SitePermutation([2,3,1])

    @test op1^-2 == SitePermutation([3,1,2])
    @test op1^-1 == SitePermutation([2,3,1])
    @test op1^0 == SitePermutation([1,2,3])
    @test op1^1 == SitePermutation([3,1,2])
    @test op1^2 == SitePermutation([2,3,1])
    @test op1^3 == SitePermutation([1,2,3])

    let
        n = -2
        @test op1^n == SitePermutation([3,1,2])
        n = -1
        @test op1^n == SitePermutation([2,3,1])
        n = 0
        @test op1^n == SitePermutation([1,2,3])
        n = 1
        @test op1^n == SitePermutation([3,1,2])
        n = 2
        @test op1^n == SitePermutation([2,3,1])
        n = 3
        @test op1^n == SitePermutation([1,2,3])
    end

    @testset "kagome" begin
        include("Kagome.jl")
        kagome = make_kagome_lattice([3 0; 0 3])

        # TranslationOperation
        t1 = TranslationOperation([1, 0])
        t1e = embed(kagome.lattice, t1)
        t2 = TranslationOperation([0, 1])
        t2e = embed(kagome.lattice, t2)

        # PointOperation
        psym = project(PointSymmetryDatabase.get(25), [1 0 0; 0 1 0])
        p6 = psym.elements[6]
        p6e = embed(kagome.lattice, p6)
        @test_throws ArgumentError embed(kagome.lattice, PointOperation([1 2; -1 1]))

        # SpaceOperation
        embed(kagome.lattice, t1 * p6)
        embed(kagome.lattice, p6 * t1)
    end

end


text:
import LinearAlgebra.schurfact!
import Base.getindex

# Schur decomposition
# redefine LinearAlgebra structure in order to support BigFloatOrComplex
struct SchurBig{Ty<:BigFloatOrComplex, S<:AbstractMatrix} <: Factorization{Ty}
    T::S
    Z::S
    values::Vector
    SchurBig{Ty,S}(T::AbstractMatrix{Ty}, Z::AbstractMatrix{Ty}, values::Vector) where {Ty,S} = new(T, Z, values)
end

function SchurBig(T::AbstractMatrix{Ty}, Z::AbstractMatrix{Ty}, values::Vector) where {Ty}
  SchurBig{Ty, typeof(T)}(T, Z, values)
end

function getindex(F::SchurBig, sym::Symbol)

  if sym == :T || sym == :Schur
    F.T
  elseif sym == :Z || sym == :vectors
    F.Z
  elseif sym == :values
    F.values
  else
    error("invalid key $(sym) for Schur")
  end
end

function schurfact!(A::StridedMatrix{T}) where {T<:BigFloatOrComplex}
  chkstride1(A)
  n = checksquare(A)
  # transform to Hessenberg form
  HF = hessfact!(A)
  Q = full(HF[:Q])
  A = HF[:H]
  separate!(A, 1, n, Q, processPart!)
  eigv = seigendiag(A, 1, size(A, 1))
  finish!(A, Q)
  SchurBig(A, Q, eigv)
end

function nonhessenberg(A::AbstractMatrix)
  n, m = size(A)
  sum = 0.0
  for i = 1:n
    for j = 1:min(i-2,m-2)
      sum += Float64(abs2(A[i,j]))
    end
  end
  for i = 2:n
    x = Float64(abs2(A[i,i-1]))
    if eltype(A) <: Complex || x < 1e-10
      # bigger elements in subdiagonal are considered block
      sum += x
    end
  end
  sqrt(sum)
end

function processPart!(A::StridedMatrix{T}, ilo::Integer, ihi::Integer, Q::AbstractMatrix) where {T<:BigFloatOrComplex}
  S = ifelse(T <: Real, Float64, Complex128)
  # 1. step: estimate eigenvalues in lower precision
  F = schurfact!(S.(view(A, ilo:ihi, ilo:ihi)))
  ev = F[:values]
  # 2. step: use low precision values ot iterate to high precision
  refineprecision!(A, ilo, ihi, Q, ev)
end

max_stars_repo_name: KlausC/LinAlgBigFloat

max_issues_repo_name: KlausC/LinAlgBigFloat

max_forks_repo_name: KlausC/LinAlgBigFloat


text:
using BackgroundSegmenter
using DataStructures

const AGE_FOR_TRANSIENCE = 2

struct Blob
    x::State
    y::State
    ratio::Float64
    area::Float64
    source::Int
    θ::Float64
    shape1::Array{Float64, 2}
    shape2::Array{Float64, 2}
    points::Array{Int, 2}
end

struct Object
    x::State
    y::State
    ratio::Float64
    area::Float64
    source::Int
    θ::Float64
    shape1::Array{Float64, 2}
    shape2::Array{Float64, 2}
    points::Array{Int, 2}
    label::Int
    age::Int
end
Object(b::Blob, label=-1, age=0) = Object(predict(b.x), predict(b.y), b.ratio, b.area, b.source, b.θ, b.shape1, b.shape2, b.points, label, age)

function update(object::Object, b::Blob, λ=0.9)
    x_velocity = smoothed_velocity_estimate(object.x, b.x, λ)
    y_velocity = smoothed_velocity_estimate(object.y, b.y, λ)
    x_observed = State(b.x.p, x_velocity)
    y_observed = State(b.y.p, y_velocity)
    Object(update(object.x, x_observed), update(object.y, y_observed),
           smooth(b.ratio, object.ratio, λ), smooth(b.area, object.area, λ),
           b.source, smooth(b.θ, object.θ, λ),
           b.shape1, b.shape2, b.points, object.label, object.age+1)
end

smooth(x, y, λ) = λ*x + (1 - λ)*y

is_transient(object::Object) = object.age < AGE_FOR_TRANSIENCE

Blob(centroid, ratio, area, source, θ) = Blob(centroid, ratio, area, source, θ,
                                              Array{Float64}(0, 0),
                                              Array{Float64}(0, 0),
                                              Array{Int}(0, 0))

function Blob(source::Int, points)
    centroid = vec(mean(points, 1))
    x = State(centroid[1], 0.0)
    y = State(centroid[2], 0.0)
    area = size(points)[1]
    if area > 1
        covariance = cov(points, 1)
    else
        covariance = [1.0 0.0; 0.0 1.0]
    end
    F = eigfact(covariance)
    ratio = (max(F[:values]...) + 0.5) / (min(F[:values]...) + 0.5)
    if F[:values][1] < F[:values][2]
        v = F[:vectors][:, 2]
        w = F[:vectors][:, 1]
    else
        v = F[:vectors][:, 1]
        w = F[:vectors][:, 2]
    end
    w .= orientation(v, w) .* w
    basis1 = hcat(v, w)
    basis2 = hcat(-v, -w)
    if v[2] == 0.0
        θ = π/2
    else
        θ = atan(v[1] / v[2])
    end
    shape1 = calc_shape(points, basis1, reshape(centroid, 1, 2))
    shape2 = calc_shape(points, basis2, reshape(centroid, 1, 2))
    return Blob(x, y, ratio, area, source, θ, shape1, shape2, points)
end

function combine(b1::Blob, b2::Blob)
    if b1.area > b2.area
        source = b1.source
    else
        source = b2.source
    end
    points = vcat(b1.points, b2.points)
    return Blob(source, points)
end

function combine(objects)
    transient = true
    source_votes = DefaultDict(Int, Float64, 0.0)
    for object in objects
        source_votes[object.source] += object.area
        if !is_transient(object)
            transient = false
        end
    end
    max_vote = 0.0
    max_source = 0
    for (source, vote) in source_votes
        if vote > max_vote
            max_vote = vote
            max_source = source
        end
    end
    points = [object.points for object in objects]
    points = vcat(points...)
    blob = Blob(max_source, points)
    return Object(blob, transient)
end

function distance(obj::Object, b::Blob)
    sqrt((obj.x.p - b.x.p)^2 + (obj.y.p - b.y.p)^2)
end

function cost(obj::Object, b::Blob)
    scale = (1.0 - pdf(obj.x, b.x.p) * pdf(obj.y, b.y.p)) * distance(obj, b)
    c = 0.0
    c += abs(obj.ratio - b.ratio) / obj.ratio
    c += abs(obj.area - b.area) / obj.area
    if obj.source != b.source
        c += 1
    end
    c += 1 - (cos(obj.θ - b.θ))
    if obj.shape1[1, :]' * b.shape1[1, :] > 0
        c += (1.0 - intersection_over_union(obj.shape1, b.shape1)) / 2
    else
        c += (1.0 - intersection_over_union(obj.shape1, b.shape2)) / 2
    end
    return c * scale
end

function orientation(x::Vector{T}, y::Vector{T}) where {T}
    x1 = zeros(3)
    y1 = zeros(3)
    x1[1:2] = x
    y1[1:2] = y
    z = cross(x1, y1)
    if z[3] > zero(T)
        return one(T)
    else
        return -one(T)
    end
end

function intersection_over_union(x, y)
    if length(x) < length(y)
        return intersection_over_union(y, x)
    end
    x = x[2:end, :]
    y = y[2:end, :]
    (a, _) = size(x)
    (b, _) = size(y)
    c1 = (a + 1.0) / 2.0
    c2 = (b + 1.0) / 2.0
    offset = convert(Int, round(c1 - c2))
    intersection = 0.0
    union = 0.0
    for j in 1:offset
        union += 2x[j, 2]
    end
    for j in (b + offset + 1):a
        union += 2x[j, 2]
    end
    for j in 1:b
        i, u = intersection_and_union(x[j+offset, :], y[j, :])
        intersection += i
        union += u
    end
    if union > 0.0
        return intersection / union
    else
        return 0.5
    end
end

function intersection_and_union(x, y)
    (c1, w1) = x
    (c2, w2) = y
    if c2 < c1
        return intersection_and_union(y, x)
    end
    union = 2.0 * (w1 + w2)
    if (c1 + w1) <= (c2 - w2)
        intersection = 0.0
    elseif (c1 - w1) >= (c2 - w2)
        intersection = 2w1
    elseif (c1 + w1) >= (c2 + w2)
        intersection = 2w2
    else
        intersection = (c1 + w1) - (c2 - w2)
    end
    union -= intersection
    return (intersection, union)
end

function composite_cost_comparison(objects, blob::Blob)
    min_cost = Inf
    best_idx = -1
    for (n, object) in enumerate(objects)
        object_cost = cost(object, blob)
        if object_cost < min_cost
            min_cost = object_cost
            best_idx = n
        end
    end
    composite = combine(objects)
    best_match = objects[best_idx]
    c = cost(composite, blob)
    if c < min_cost
        return -1
    else
        return best_idx
    end
end

function calc_shape(points, basis, centroid)
    deskewed_points = (points .- centroid) * basis
    min_idx = Inf
    max_idx = -Inf
    (n, m) = size(deskewed_points)
    ranges = Dict{Int, Vector{Float64}}()
    for i in 1:n
        idx = convert(Int, round(deskewed_points[i, 1]))
        if idx < min_idx
            min_idx = idx
        end
        if idx > max_idx
            max_idx = idx
        end
        if !haskey(ranges, idx)
            ranges[idx] = [deskewed_points[i, 2]]
        else
            push!(ranges[idx], deskewed_points[i, 2])
        end
    end
    shape = Array{Float64}(max_idx - min_idx + 2, 2)
    shape[1, :] = basis[:, 1]'
    for (k, idx) in enumerate(min_idx:max_idx)
        if haskey(ranges, idx)
            b = max(ranges[idx]...)
            a = min(ranges[idx]...)
            shape[k+1, 1] = (b + a) / 2
            shape[k+1, 2] = (b - a) / 2
        else
            shape[k+1, 1] = 0.0
            shape[k+1, 2] = 0.0
        end
    end
    return shape
end

function form_blobs(w::Array{T, 3}) where {T}
    (n, m, t) = size(w)
    blobs = Vector{Vector{Blob}}(t)
    for j in 1:t
        blobs[j] = form_blobs(view(w, :, :, j))
    end
    return blobs
end

function form_blobs(w::AbstractArray{T, 2}) where {T}
    blobs = Vector{Blob}(0)
    C = components(w)
    for (label, points) in C
        push!(blobs, Blob(label, points))
    end
    return blobs
end


text:
using PEGParser

@grammar lisp_grammar begin
  float_with_dot_matcher = r"[-+]?[0-9]*\.[0-9]+([eE][-+]?[0-9]+)?[Ff]" { float32(_0[1:end-1]) }
  float_no_dot_matcher   = r"[-+]?[0-9]*[0-9]+([eE][-+]?[0-9]+)?[Ff]" { float32(_0[1:end-1]) }
  white_space            = -r"([\s\n\r]*(?<!\\);[^\n\r$]+[\n\r\s$]*|[\s\n\r]+)"

  doubley  = r"[-+]?[0-9]*\.?[0-9]+([eE][-+]?[0-9]+)?[dD]" { float64(_0[1:end-1]) }
  floaty   = (float_with_dot_matcher | float_no_dot_matcher) { _1 }
  inty     = r"[-+]?\d+" { int(_0) }
  uchary   = r"\\(u[\da-fA-F]{4})" { begin x = unescape_string(_0); x[chr2ind(x, 1)] end }
  achary   = r"\\[0-7]{3}" { unescape_string(_0)[1] }
  chary    = r"\\." { _0[2] }
  stringy  = r"(?<!\\)\".*?(?<!\\)\"" { _0[2:end-1] } #_0[2:end-1] }
  booly    = r"(true|false)"i { _0 == "true" ? true : false }
  symboly  = r"[^\d(){}#'`,@~;~\[\]^\s][^\s()#'`,@~;^{}~\[\]]*" { symbol(_0) }
  sexpr0   = ("(" + ?(white_space) + ")") { {{}} }
  sexpr1   = ("(" + list(expr, white_space) + ")") { { vcat(_2.children...) } }
  sexpr    = (sexpr1 | sexpr0) { children }
  hashy    = ("#{" + list(expr, white_space) + "}") { Set(vcat(_2.children...)...) }
  curly    = ("{" + list(expr, white_space) + "}") { [ _2.children[i][1] => _2.children[i+1][1] for i = 1:2:length(_2.children) ] }
  bracket1 = ("[" + list(expr, white_space) + "]") { { vcat(_2.children...) } }
  bracket0 = ("[" + ?(white_space) + "]") { {{}} }
  bracket  = (bracket0 | bracket1) { children }
  quot     = ("'" + expr) { {[ :quote _2[1] ]} }
  quasi    = ("`" + expr) { {[ :quasi _2 ]} }
  tildeseq = ("~@" + expr) { {[ :splice_seq, { _2[1] } ]} }
  tilde    = ("~" + expr) { {{ :splice, _2[1] }} }
  # expr    = (double | float | int | uchar | achar | char | string | bools | symbol | sexpr |
  #            hash | quot | curly | bracket | quasi | tilde)
  # start   = expr
  expr    = (doubley | floaty | inty | uchary | achary | chary | stringy | booly | symboly | sexpr | hashy | curly | bracket | quot | quasi | tildeseq | tilde) { children }
  start   = expr { _1 }
end

function read(input)
  ast, pos, err = parse(lisp_grammar, input)
  if ast == nothing
    println("pos: " * string(pos))
    println("err: " * string(err))
    return nothing
  else
    #clean = transform(grammar_transform, ast)#, ignore = [:white_space, :float_with_dot_matcher, :float_no_dot_matcher])
    #return clean
    return ast
  end
end


text:
function hello_to(name::AbstractString)
    num_greetings = length(GREETINGS)
    rand_greeting = round(Int, (num_greetings-1)*rand()) + 1
    string(GREETINGS[rand_greeting], " ", name)
end


text:
@testset "parameters" begin
    function objective_function(unflatten, flat_θ::Vector{<:Real})
        θ = value(unflatten(flat_θ))
        return abs2(θ.a) + abs2(θ.b)
    end

    # This is more of a worked example.
    @testset "Integration" begin
        θ0 = (a=5.0, b=4.0)
        flat_parameters, unflatten = flatten(θ0)

        results = Optim.optimize(
            θ -> objective_function(unflatten, θ),
            θ -> only(Zygote.gradient(θ -> objective_function(unflatten, θ), θ)),
            flat_parameters,
            BFGS(),
            Optim.Options();
            inplace=false,
        )

        # Check that it's successfully optimised.
        @test objective_function(unflatten, results.minimizer) < 1e-12
    end

    @testset "Other Integration" begin
        θ0 = (a=5.0, b=fixed(4.0))
        flat_parameters, unflatten = flatten(θ0)

        results = Optim.optimize(
            θ -> objective_function(unflatten, θ),
            θ -> only(Zygote.gradient(θ -> objective_function(unflatten, θ), θ)),
            flat_parameters,
            BFGS(),
            Optim.Options();
            inplace=false,
        )

        # Check that it's successfully optimised.
        @test value(unflatten(results.minimizer).b) == 4.0
    end

    function other_objective_function(unflatten, flat_θ::Vector{<:Real})
        X = value(unflatten(flat_θ))
        return -logpdf(X, -1.0) - logpdf(X, 1.0)
    end

    @testset "Normal" begin
        θ0 = deferred(Normal, randn(), positive(1.0))
        flat_parameters, unflatten = flatten(θ0)

        results = Optim.optimize(
            θ -> other_objective_function(unflatten, θ),
            θ -> only(Zygote.gradient(θ -> other_objective_function(unflatten, θ), θ)),
            flat_parameters,
            BFGS(),
            Optim.Options();
            inplace=false,
        )

        # Check that it's successfully optimised.
        @test mean(value(unflatten(results.minimizer))) ≈ 0 atol = 1e-7
    end

    @testset "value_flatten" begin
        x = (ones(3), fixed(5.0), (a=fixed(5.0), b=[6.0, 2.1]))
        v, unflatten = value_flatten(x)

        @test length(v) == 5
        @test unflatten(v) == ParameterHandling.value(x)
    end
end


text:
# ------------------------------------------------------------------
# Licensed under the MIT License. See LICENSE in the project root.
# ------------------------------------------------------------------

"""
    NuggetEffect(nugget=n)

A pure nugget effect variogram with nugget `n`.
"""
struct NuggetEffect{V} <: Variogram
  nugget::V
end

NuggetEffect(; nugget=1.0) = NuggetEffect(nugget)

(γ::NuggetEffect)(h) = (h > 0) * γ.nugget

function (γ::NuggetEffect)(u::Point, v::Point)
  d = Euclidean()
  x = coordinates(u)
  y = coordinates(v)
  h = evaluate(d, x, y)
  γ(h)
end

Base.range(::NuggetEffect{T}) where {T} = zero(T)

sill(γ::NuggetEffect) = γ.nugget

isstationary(::Type{<:NuggetEffect}) = true

isisotropic(::NuggetEffect) = true


text:
using PyCall
py"print('Hello World from Python')"
py"2+2"
A = 15
py"print"("A = $A")

os = pyimport("os")

cwd = os.getcwd()
py"print"("cwd = $cwd")
os.path.basename(cwd)
os.path.dirname(cwd)


text:
using OrdinaryDiffEq, Distributions,
      DiffEqUncertainty, Test, Quadrature, Cubature

function f(du,u,p,t)
  @inbounds begin
    du[1] = dx = p[1]*u[1] - u[1]*u[2]
    du[2] = dy = -3*u[2] + u[1]*u[2]
  end
  nothing
end

u0 = [1.0;1.0]
tspan = (0.0,10.0)
p = [1.5]
prob = ODEProblem(f,u0,tspan,p)
sol = solve(remake(prob,u0=u0),Tsit5())
cost(sol) = sum(max(x[1]-6,0) for x in sol.u)
u0s = [Uniform(0.25,5.5),Uniform(0.25,5.5)]
ps  = [Uniform(0.5,2.0)]
@time sol = koopman_expectation(cost,u0s,ps,prob,Tsit5();iabstol=1e-3,ireltol=1e-3,maxiters=1000,saveat=0.1)
c1, e1 = sol.u, sol.resid
@time sol = koopman_expectation(cost,u0s,ps,prob,Tsit5(),EnsembleThreads();quadalg=CubatureJLh(),
                         batch=1000,iabstol=1e-3,ireltol=1e-3,
                         maxiters=2000,saveat=0.1)
c2, e2 = sol.u, sol.resid
@test abs(c1 - c2) < 0.1

@time c3 = montecarlo_expectation(cost,u0s,ps,prob,Tsit5(),EnsembleThreads();trajectories=100000,saveat=0.1)
@test abs(c1 - c3) < 0.1

#=
using DiffEqGPU
@time sol = koopman_expectation(cost,u0s,ps,prob,Tsit5(),EnsembleGPUArray();quadalg=CubatureJLh(),
                         batch=1000,iabstol=1e-3,ireltol=1e-3,
                         maxiters=2000,saveat=0.1)
c2, e2 = sol.u, sol.resid
@test abs(c1 - c2) < 0.1

@time c3 = montecarlo_expectation(cost,u0s,ps,prob,Tsit5(),EnsembleGPUArray();trajectories=100000,saveat=0.1)
@test abs(c1 - c3) < 0.1
=#


text:
#! format: off

using SIIPExamples
pkgpath = pkgdir(SIIPExamples)
using Dates

using PowerSystems
using PowerSimulations
using PowerSystemCaseBuilder

using HiGHS # mip solver
solver = optimizer_with_attributes(HiGHS.Optimizer, "mip_rel_gap" => 0.5)

sys = build_system(PSYTestSystems, "tamu_ACTIVSg2000_sys")
transform_single_time_series!(sys, 24, Dates.Hour(24));

sim_folder = mkpath(joinpath(pkgpath, "TAMU-sim"))
models = SimulationModels(
    decision_models = [
        DecisionModel(template_unit_commitment(), sys, name = "UC", optimizer = solver),
    ],
)

sim = Simulation(
    name = "TAMU-test",
    steps = 3,
    models = models,
    sequence = SimulationSequence(
        models = models,
        ini_cond_chronology = InterProblemChronology(),
    ),
    simulation_folder = sim_folder,
)

build!(sim)

execute!(sim)

sim_results = SimulationResults(sim);
uc_results = get_decision_problem_results(sim_results, "UC")
read_realized_variable(uc_results, "OnVariable__ThermalStandard")

# This file was generated using Literate.jl, https://github.com/fredrikekre/Literate.jl

max_stars_repo_name: NREL-SIIP/Examples

max_issues_repo_name: NREL-SIIP/Examples

max_forks_repo_name: NREL-SIIP/Examples


text:
using MechanicalSketch
import MechanicalSketch: sethue, background, foil_spline_local, boxtopleft, boxtopright
import MechanicalSketch: BoundingBox, m, mm, °, polyrotate!, polymove!, color_with_lumin
import MechanicalSketch: O, EM, poly, lumin, get_current_lumin, get_current_RGB, gsave, grestore
import MechanicalSketch: empty_figure, PALETTE,  dimension_aligned, color_from_palette, Point

let
BACKCOLOR = color_with_lumin(PALETTE[8], 70)
function restart()
    empty_figure(filename = joinpath(@__DIR__, "test_6.png"))
    background(BACKCOLOR)
    sethue(PALETTE[7])
end
restart()
O + (0.0m, 1.0m)
SC = 15
fsp = foil_spline_local(l = SC*1000mm, t = SC * 58.4mm, c= SC * 30.76mm);
# Draw the thing, with outline

poly(fsp, :fill)
sethue(color_from_palette("indigo"))
poly(fsp, :stroke)
# Pretend we don't know the size, measure and dimension
b = BoundingBox(fsp)
dimension_aligned(boxtopleft(b), boxtopright(b))
# Test to "measure" with more precision than what we get with boundingbox
p1 = Point(-0.3 *SC * 1000mm, 0.0mm)
p2 = Point(0.7 *SC * 1000mm, 0.0mm)
dimension_aligned(p1, p2; unit = mm, offset = -8EM, digits = 0)

gsave()
for i = 0:9
    sethue(PALETTE[i + 1])
    f = SC .* foil_spline_local();
    polyrotate!(f, i * °)
    polymove!(f,  O, O + (-5m, -2m))
    poly(f, :fill)
    bckco = BACKCOLOR
    luminback = lumin(bckco)
    luminfront = get_current_lumin()
    deltalumin = luminback - luminfront
    contrastcol = color_with_lumin(get_current_RGB(), luminfront - deltalumin)
    sethue(contrastcol)
    poly(f, :stroke)
end
grestore()

MechanicalSketch.finish()
nothing
end


text:
module MultiNonLinearFEOperatorsTests

using Test
using Gridap
import Gridap: ∇

# Define manufactured functions
u1fun(x) = x[1] + x[2]
u2fun(x) = x[1] - x[2]

u1fun_grad(x) = VectorValue(1.0,1.0)
u2fun_grad(x) = VectorValue(1.0,-1.0)

∇(::typeof(u1fun)) = u1fun_grad
∇(::typeof(u2fun)) = u2fun_grad

b1fun(x) = u2fun(x) -(3.0*x[1]+x[2]+1.0)
b2fun(x) = 0.0

@law ν(x,u1) = (u1+1.0)*x[1]
@law dν(x,du1) = du1*x[1]

# Construct the discrete model
model = CartesianDiscreteModel(domain=(0.0,1.0,0.0,1.0), partition=(4,4))

# Construct the FEspace
order = 2
diritag = "boundary"
fespace = H1ConformingFESpace(Float64,model,order,diritag)

# Define test and trial
V1 = TestFESpace(fespace)
V2 = V1
V = [V1, V2]

U1 = TrialFESpace(fespace,u1fun)
U2 = TrialFESpace(fespace,u2fun)
U = [U1, U2]

# Define integration mesh and quadrature
trian = Triangulation(model)
quad = CellQuadrature(trian,degree=3*order-1)

# Define cell field describing the source term
b1field = CellField(trian,b1fun)
b2field = CellField(trian,b2fun)

# Define residual and jacobian
a(u,v,du) = inner(∇(v[1]),ν(u[1])*∇(du[1])) + inner(v[1],du[2]) + inner(∇(v[2]),∇(du[2]))
b(v) = inner(v[1],b1field) + inner(v[2],b2field)

res(u,v) = a(u,v,u) - b(v)
jac(u,v,du) = a(u,v,du) + inner(∇(v[1]),dν(du[1])*∇(u[1]))

# Define Assembler
assem = SparseMatrixAssembler(V,U)

# Define the FEOperator
op = NonLinearFEOperator(res,jac,V,U,assem,trian,quad)

# Define the FESolver
ls = LUSolver()
tol = 1.e-10
maxiters = 20
nls = NewtonRaphsonSolver(ls,tol,maxiters)
solver = NonLinearFESolver(nls)

# Solve!
uh = solve(solver,op)

# Define exact solution and error
u1 = CellField(trian,u1fun)
e1 = u1 - uh[1]

u2 = CellField(trian,u2fun)
e2 = u2 - uh[2]

# Define norms to measure the error
l2(u) = inner(u,u)
h1(u) = inner(∇(u),∇(u)) + l2(u)

# Compute errors
e1l2 = sqrt(sum( integrate(l2(e1),trian,quad) ))
e1h1 = sqrt(sum( integrate(h1(e1),trian,quad) ))

e2l2 = sqrt(sum( integrate(l2(e2),trian,quad) ))
e2h1 = sqrt(sum( integrate(h1(e2),trian,quad) ))

@test e1l2 < 1.e-8
@test e1h1 < 1.e-8

@test e2l2 < 1.e-8
@test e2h1 < 1.e-8

# Further tests

t_Ω = NonLinearFETerm(res,jac,trian,quad)
op = NonLinearFEOperator(V,U,t_Ω)

uh = solve(solver,op)

end