text:
#
# Single particle von Neumann entanglement entropy
#

function entropy(n::Number)
  (n ≤ 0 || n ≥ 1) && return 0
  return -(n*log(n) + (1-n)*log(1-n))
end

entropy(ns::Vector{Float64}) = sum(entropy, ns)

#
# Linear Algebra tools
#

"""
    frobenius_distance(M1::AbstractMatrix, M2::AbstractMatrix)

Computes the Frobenius distance `√tr((M1-M2)'*(M1-M2))`.
"""
frobenius_distance(M1::AbstractMatrix, M2::AbstractMatrix) =
  sqrt(abs(tr(M1'M1) + tr(M2'M2) - tr(M1'M2) - tr(M2'M1)))

#
# Rotations
#

struct Circuit{T} <: LinearAlgebra.AbstractRotation{T}
  rotations::Vector{Givens{T}}
end

Base.adjoint(R::Circuit) = Adjoint(R)

Base.copy(aR::Adjoint{<:Any,Circuit{T}}) where {T} = Circuit{T}(reverse!([r' for r in aR.parent.rotations]))

function LinearAlgebra.lmul!(G::Givens, R::Circuit)
  push!(R.rotations, G)
  return R
end

function LinearAlgebra.lmul!(R::Circuit, A::AbstractArray)
  @inbounds for i in 1:length(R.rotations)
      lmul!(R.rotations[i], A)
  end
  return A
end

function LinearAlgebra.rmul!(A::AbstractMatrix, adjR::Adjoint{<:Any,<:Circuit})
  R = adjR.parent
  @inbounds for i in 1:length(R.rotations)
      rmul!(A, adjoint(R.rotations[i]))
  end
  return A
end

Base.:*(g1::Circuit, g2::Circuit) = Circuit(vcat(g2.rotations, g1.rotations))
LinearAlgebra.lmul!(g1::Circuit, g2::Circuit) = append!(g2.rotations, g1.rotations)

Base.:*(A::Circuit, B::Union{ <: Hermitian, <: Diagonal}) = A * convert(Matrix, B)
Base.:*(A::Adjoint{ <: Any, <: Circuit}, B::Hermitian) = copy(A) * convert(Matrix, B)
Base.:*(A::Adjoint{ <: Any, <: Circuit}, B::Diagonal) = copy(A) * convert(Matrix, B)
Base.:*(A::Adjoint{<:Any, <: AbstractVector}, B::Adjoint{ <:Any, <: Circuit}) = convert(Matrix, A) * B

function LinearAlgebra.rmul!(A::AbstractMatrix, R::Circuit)
  @inbounds for i in reverse(1:length(R.rotations))
      rmul!(A, R.rotations[i])
  end
  return A
end

function shift!(G::Circuit, i::Int)
  for (n, g) in enumerate(G.rotations)
    G.rotations[n] = Givens(g.i1+i, g.i2+i, g.c, g.s)
  end
  return G
end

ngates(G::Circuit) = length(G.rotations)

#
# Free fermion tools
#

# Make a hopping Hamiltonian from quadratic Hamiltonian
function hopping_hamiltonian(ampo::AutoMPO)
  nterms = length(ampo.data)
  coefs = Vector{Number}(undef, nterms)
  sites = Vector{Tuple{Int, Int}}(undef, nterms)
  nsites = 0
  for n in 1:nterms
    term = ampo.data[n]
    coef = isreal(term.coef) ? real(term.coef) : term.coef
    coefs[n] = coef
    ops = term.ops
    length(ops) != 2 && error("Must create hopping Hamiltonian from quadratic Hamiltonian")
    sites[n] = ntuple(n -> ops[n].site, Val(2))
    nsites = max(nsites, maximum(sites[n]))
  end
  ElT = all(isreal(coefs)) ? Float64 : ComplexF64
  h = zeros(ElT, nsites, nsites)
  for n in 1:nterms
    h[sites[n]...] = coefs[n]
  end
  return h
end

# Make a combined hopping Hamiltonian for spin up and down
function hopping_hamiltonian(ampo_up::AutoMPO, ampo_dn::AutoMPO)
  h_up = hopping_hamiltonian(ampo_up)
  h_dn = hopping_hamiltonian(ampo_dn)
  @assert size(h_up) == size(h_dn)
  N = size(h_up, 1)
  ElT = promote_type(eltype(h_up), eltype(h_dn))
  h = zeros(ElT, 2*N, 2*N)
  for i in 1:2*N, j in 1:2*N
    if isodd(i) && isodd(j)
      i_up, j_up = (i+1)÷2, (j+1)÷2
      h[i, j] = h_up[i_up, j_up]
    elseif iseven(i) && iseven(j)
      i_dn, j_dn = i÷2, j÷2
      h[i, j] = h_dn[i_dn, j_dn]
    end
  end
  return Hermitian(h)
end

# Make a Slater determinant matrix from a hopping Hamiltonian
# h with Nf fermions.
function slater_determinant_matrix(h::AbstractMatrix, Nf::Int)
  _, u = eigen(h)
  return u[:, 1:Nf]
end

function Base.:+(a1::AutoMPO, a2::AutoMPO)
  return AutoMPO(vcat(a1.data, a2.data))
end

#
# Correlation matrix diagonalization
#

"""
    givens_rotations(v::AbstractVector)

For a vector `v`, return the `length(v)-1`
Givens rotations `g` and the norm `r` such that:
```julia
g * v ≈ r * [n == 1 ? 1 : 0 for n in 1:length(v)]
```
"""
function givens_rotations(v::AbstractVector{ElT}) where {ElT}
  N = length(v)
  gs = Circuit{ElT}([])
  r = v[1]
  for n in reverse(1:N-1)
    g, r = givens(v, n, n+1)
    v = g * v
    lmul!(g, gs)
  end
  return gs, r
end

"""
    correlation_matrix_to_gmps(Λ::AbstractMatrix{ElT}; eigval_cutoff::Float64 = 1e-8, maxblocksize::Int = size(Λ0, 1))

Diagonalize a correlation matrix, returning the eigenvalues and eigenvectors
stored in a structure as a set of Givens rotations.

The correlation matrix should be Hermitian, and will be treated as if it itensor
in the algorithm.
"""
function correlation_matrix_to_gmps(Λ0::AbstractMatrix{ElT}; eigval_cutoff::Float64 = 1e-8,
                                                             maxblocksize::Int = size(Λ0, 1)) where {ElT <: Number}
  Λ = Hermitian(Λ0)
  N = size(Λ, 1)
  V = Circuit{ElT}([])
  ns = Vector{real(ElT)}(undef, N)
  err_tot = 0.0
  for i in 1:N
    blocksize = 0
    n = 0.0
    err = 0.0
    p = Int[]
    uB = 0.0
    for blocksize in 1:maxblocksize
      j = min(i + blocksize, N)
      ΛB = @view Λ[i:j, i:j]
      nB, uB = eigen(Hermitian(ΛB))
      p = sortperm(nB; by = entropy)
      n = nB[p[1]]
      err = min(n, 1-n)
      err ≤ eigval_cutoff && break
    end
    err_tot += err
    ns[i] = n
    v = @view uB[:, p[1]]
    g, _ = givens_rotations(v)
    shift!(g, i-1)
    # In-place version of:
    # V = g * V
    lmul!(g, V)
    Λ = Hermitian(g * Λ * g')
  end
  return ns, V
end

function slater_determinant_to_gmps(Φ::AbstractMatrix; kwargs...)
  return correlation_matrix_to_gmps(conj(Φ) * transpose(Φ); kwargs...)
end

#
# Turn circuit into MPS
#

function ITensors.ITensor(u::Givens, s1::Index, s2::Index)
  U = [1    0   0 0
       0  u.c u.s 0
       0 -conj(u.s) u.c 0
       0    0   0 1]
  return itensor(U, s2', s1', dag(s2), dag(s1))
end

function ITensors.ITensor(sites::Vector{<:Index}, u::Givens)
  s1 = sites[u.i1]
  s2 = sites[u.i2]
  return ITensor(u, s1, s2)
end

"""
    MPS(sites::Vector{<:Index}, state, U::Vector{<:ITensor}; kwargs...)

Return an MPS with site indices `sites` by applying the circuit `U` to the starting state `state`.
"""
function ITensors.MPS(sites::Vector{<:Index}, state, U::Vector{<:ITensor}; kwargs...)
  return apply(U, productMPS(sites, state); kwargs...)
end

function isspinful(s::Index)
  !hasqns(s) && return false
  return all(qnblock -> ITensors.hasname(qn(qnblock), ITensors.QNVal("Sz",0)), space(s))
end

function isspinful(s::Vector{ <: Index})
  return all(isspinful, s)
end

"""
    correlation_matrix_to_mps(Λ::AbstractMatrix{ElT}; eigval_cutoff::Float64 = 1e-8,
                                                      maxblocksize::Int = size(Λ, 1),
                                                      kwargs...)

Return an approximation to the state represented by the correlation matrix as
a matrix product state (MPS).

The correlation matrix should correspond to a pure state (have all eigenvalues
of zero or one).
"""
function correlation_matrix_to_mps(s::Vector{<:Index}, Λ::AbstractMatrix;
                                   eigval_cutoff::Float64 = 1e-8, maxblocksize::Int = size(Λ, 1), kwargs...)
  @assert size(Λ, 1) == size(Λ, 2)
  ns, C = correlation_matrix_to_gmps(Λ; eigval_cutoff = eigval_cutoff, maxblocksize = maxblocksize)
  if all(hastags("Fermion"), s)
    U = [ITensor(s, g) for g in reverse(C.rotations)]
    ψ = MPS(s, n -> round(Int, ns[n]) + 1, U; kwargs...)
  elseif all(hastags("Electron"), s)
    isodd(length(s)) && error("For Electron type, must have even number of sites of alternating up and down spins.")
    N = length(s)
    if isspinful(s)
      error("correlation_matrix_to_mps(Λ::AbstractMatrix) currently only supports spinless Fermions or Electrons that do not conserve Sz. Use correlation_matrix_to_mps(Λ_up::AbstractMatrix, Λ_dn::AbstractMatrix) to use spinful Fermions/Electrons.")
    else
      sf = siteinds("Fermion", 2*N; conserve_qns = true)
    end
    U = [ITensor(sf, g) for g in reverse(C.rotations)]
    ψf = MPS(sf, n -> round(Int, ns[n]) + 1, U; kwargs...)
    ψ = MPS(N)
    for n in 1:N
      i, j = 2*n-1, 2*n
      C = combiner(sf[i], sf[j])
      c = combinedind(C)
      ψ[n] = ψf[i] * ψf[j] * C
      ψ[n] *= δ(dag(c), s[n])
    end
  else
    error("All sites must be Fermion or Electron type.")
  end
  return ψ
end

function slater_determinant_to_mps(s::Vector{<:Index}, Φ::AbstractMatrix; kwargs...)
  return correlation_matrix_to_mps(s, conj(Φ) * transpose(Φ); kwargs...)
end

function slater_determinant_to_mps(s::Vector{<:Index}, Φ_up::AbstractMatrix, Φ_dn::AbstractMatrix; kwargs...)
  return correlation_matrix_to_mps(s, conj(Φ_up) * transpose(Φ_up), conj(Φ_dn) * transpose(Φ_dn); kwargs...)
end

function mapindex(f::Function, C::Circuit)
  return Circuit(mapindex.(f, C.rotations))
end

function mapindex(f::Function, g::Givens)
  return Givens(f(g.i1), f(g.i2), g.c, g.s)
end

function replace_index_itensor(i1::Index, i2::Index)
  # Should be diagonal within each block
  # TODO: make a special constructor ITensor(1.0, ::Index, ::Index) (this currently doesn't work
  # for QN ITensors)
  # TODO: make a special constructor diagBlockITensor?
  t = ITensor(dag(i1), i2)
  # TODO: use nnz
  ITensors.data(t) .= 1
  return t
end

function correlation_matrix_to_mps(s::Vector{<:Index}, Λ_up::AbstractMatrix, Λ_dn::AbstractMatrix;
                                   eigval_cutoff::Float64 = 1e-8,
                                   maxblocksize::Int = min(size(Λ_up, 1), size(Λ_dn, 1)),
                                   kwargs...)
  @assert size(Λ_up, 1) == size(Λ_up, 2)
  @assert size(Λ_dn, 1) == size(Λ_dn, 2)
  N_up = size(Λ_up, 1)
  N_dn = size(Λ_dn, 1)
  @assert N_up == N_dn
  # Total number of fermion sites
  N = N_up + N_dn
  ns_up, C_up = correlation_matrix_to_gmps(Λ_up; eigval_cutoff = eigval_cutoff, maxblocksize = maxblocksize)
  ns_dn, C_dn = correlation_matrix_to_gmps(Λ_dn; eigval_cutoff = eigval_cutoff, maxblocksize = maxblocksize)
  # map the up spins to the odd sites and the even spins to the even sites
  C_up = mapindex(n -> 2n-1, C_up)
  C_dn = mapindex(n -> 2n, C_dn)
  C = Circuit(collect(Iterators.flatten(zip(C_up.rotations, C_dn.rotations))))
  ns = collect(Iterators.flatten(zip(ns_up, ns_dn)))
  if all(hastags("Fermion"), s)
    @assert length(s) == N
    U = [ITensor(s, g) for g in reverse(C.rotations)]
    ψ = MPS(s, n -> round(Int, ns[n]) + 1, U; kwargs...)
  elseif all(hastags("Electron"), s)
    @assert length(s) == N_up
    @assert length(s) == N_dn
    if isspinful(s)
      space_up = [QN(("Nf", 0, -1), ("Sz", 0)) => 1,
                  QN(("Nf", 1, -1), ("Sz", 1)) => 1]
      space_dn = [QN(("Nf", 0, -1), ("Sz", 0)) => 1,
                  QN(("Nf", 1, -1), ("Sz", -1)) => 1]
      sf_up = [Index(space_up, "Fermion,Site,n=$(2n-1)") for n in 1:N_up]
      sf_dn = [Index(space_dn, "Fermion,Site,n=$(2n)") for n in 1:N_dn]
      sf = collect(Iterators.flatten(zip(sf_up, sf_dn)))
    else
      sf = siteinds("Fermion", N; conserve_qns = true, conserve_sz = false)
    end
    U = [ITensor(sf, g) for g in reverse(C.rotations)]
    ψf = MPS(sf, n -> round(Int, ns[n]) + 1, U; kwargs...)
    ψ = MPS(N_up)
    for n in 1:N_up
      i, j = 2*n-1, 2*n
      C = combiner(sf[i], sf[j])
      c = combinedind(C)
      ψ[n] = ψf[i] * ψf[j] * C
      ψ[n] *= replace_index_itensor(c, s[n])
    end
  else
    error("All sites must be Fermion or Electron type.")
  end
  return ψ
end


text:
using Documenter, Embeddings

makedocs(;
    modules=[Embeddings],
    format=:html,
    pages=[
        "Home" => "index.md",
    ],
    repo="https://github.com/JuliaText/Embeddings.jl/blob/{commit}{path}#L{line}",
    sitename="Embeddings.jl",
    authors="Lyndon White (aka oxinabox)",
    assets=[],
)

deploydocs(;
    repo="github.com/JuliaText/Embeddings.jl",
    target="build",
    julia="0.6",
    deps=nothing,
    make=nothing,
)

ext: jl


text:
using Test
using LinearAlgebra

import Pkg
Pkg.add("FFTW")
using FFTW

Fₘ = fft(I(4), 1)

imag(Fₘ)

function CreateLM(N=2)
  L = 2 * ones(1, N) - 2 * cos.(1:N)' * π / (N + 1)
  LM = ones(N, 1) * L + L' * ones(1, N)
end

CreateLM(2)

# Packing 2D Vector into Matrices
N = 4
L = 2 * ones(1, N) - 2 * cos.(1:N)' * π / (N + 1)
LM = ones(N, 1) * L + L' * ones(1, N)


text:
struct MyPoly
    c
    var
end

p = MyPoly([3,-5,2], 'x')

struct MyPoly
    c
    var
    MyPoly(c, var='x') = new(c, var)
end

p.c

p.var

resize!(p.c, 4)
p.c[1:4] = [-1,3,-5,2]
p

p = MyPoly([-1,3,-5,2], p.var)

function times_x(p)
   return MyPoly(vcat(p.c,0), p.var)
end

times_x(p)

function degree(p::MyPoly)
    ix1 = findfirst(p.c .!= 0)
    if ix1 == nothing
        return 0
    else
        return length(p.c) - ix1
    end
end

println(degree(MyPoly([0,0,0,0,0])))
println(degree(MyPoly([0,0,0,0,1])))
println(degree(MyPoly([0,0,0,1,0])))
println(degree(MyPoly([1,0,0,0,0])))

function degree(p::Int)
    println("degree function called with Int argument")
end

function degree(p)
    println("degree function called with any other argument")
end

degree(1.234)
degree([1,2])
degree(-123)
degree(MyPoly([1,2,3]))

function Base.show(io::IO, p::MyPoly)
    d = degree(p)
    print(io, "MyPoly(")
    for k = d:-1:0
        coeff = p.c[end-k]
        if coeff == 0 && d > 0
            continue
        end
        if k < d
            if isa(coeff, Real)
                if coeff > 0
                    print(io, " + ")
                else
                    print(io, " - ")
                end
                coeff = abs(coeff)
            else
                print(io, " + ")
            end
        end
        if isa(coeff, Real)
            print(io, coeff)
        else
            print(io, "($coeff)")
        end
        if k == 0
            continue
        end
        print(io, "⋅", p.var)
        if k > 1
            print(io, "^", k)
        end
    end
    print(io, ")")
end

p

MyPoly([-1.234,0,0,0,4.321], 's')

function (p::MyPoly)(x)
    d = degree(p)
    v = p.c[end-d]
    for cc = p.c[end-d+1:end]
        v = v*x + cc
    end
    return v
end

println(p(1))
println(p(1.234))
println(p.([1,-2,1.3]))  # Note: Broadcasting automatically defined

using PyPlot
function PyPlot.plot(p::MyPoly, xlim=[-2,2])
    xx = collect(range(xlim[1], xlim[2], length=100))
    plot(xx, p.(xx))
    xlabel(string(p.var))
end

p = MyPoly([1,1,-5,-5,4,4])
plot(p)
p

function Base.:+(p1::MyPoly, p2::MyPoly)
    if p1.var != p2.var
        error("Cannot add polynomials of different independent variables.")
    end
    d1 = degree(p1)
    d2 = degree(p2)
    d = max(d1,d2)
    c = [fill(0, d-d1); p1.c] + [fill(0, d-d2); p2.c]
    return MyPoly(c, p1.var)
end

println(p)
println(p + p)
println(p + MyPoly([1.1,2.2]))
println(p + MyPoly([1], 's')) # will trigger our error message

function Base.:-(p1::MyPoly, p2::MyPoly)
    return p1 + MyPoly(-p2.c)
end

function Base.:*(a, p::MyPoly)
    newc = a * p.c
    return MyPoly(newc, p.var)
end

function Base.:*(p::MyPoly, a)
    return a*p
end

p1 = 0.4p
p2 = p1 - .3*MyPoly([-2,1,1])
p3 = -1p2 + p
plot.([p1,p2,p3]);

p = MyPoly([1, -2//3, 6//7])

p.([1, -7//2])

p = MyPoly([1, im, -1, -im, 1])

p.([0, im, -1 - im])

p = MyPoly(collect(-1.5:3.5))

p(BigFloat(-π))


text:
"""
A small 2x2 GridNetwork extended in time by 2 steps
All edges, including temporal ones, have unit costs
"""

using QuNet

g = GridNetwork(2,2)
smalltemp = QuNet.TemporalGraph(g, 2, memory_costs=unit_costvector())

max_stars_repo_name: FalafelGood/QuNet

max_issues_repo_name: FalafelGood/QuNet

max_forks_repo_name: FalafelGood/QuNet


text:
using Compat
using Knet
using Base.Test
using Base.Test: Success, Failure, Error
import Base.Test: default_handler
include("isapprox.jl")
if Knet.GPU
    eval(Expr(:using,:CUDArt))
    eval(Expr(:using,:CUDArt,:ContiguousArray))
else
    typealias ContiguousArray{T} Array{T}
end

gnet0 = net0 = x0 = z0 = nothing

# Uncomment this if you want lots of messages:
# default_handler(r::Success) = info("$(r.expr)")
# default_handler(r::Failure) = warn("FAIL: $(r.expr)")
# default_handler(r::Error)   = warn("$(r.err): $(r.expr)")

function gradcheck(net, x1, z1, w, dw; iter=10,
                   epsilon=eps(eltype(w))^(1/4),
                   delta=(eltype(w)==Float64) ? 1e-4 : 1e-3)
    for i=1:iter
        r = (iter==length(w) ? i : rand(1:length(w)))
        wr0 = w[r]
        wr1 = wr0 - delta
        wr2 = wr0 + delta
        # Do not cross 0 for softloss
        (wr0>=0) && (wr1<0) && (wr1=0)
        (wr0<0) && (wr2>0) && (wr2=0)
        w[r] = wr1; loss1 = getloss(net, x1, z1)
        w[r] = wr2; loss2 = getloss(net, x1, z1)
        w[r] = wr0
        dwr = (loss2 - loss1) / (wr2 - wr1)
        # @show (dw[r], dwr)
        if !isapprox(dw[r], dwr; rtol=0.1)
            @show (:gradcheck, dw[r], dwr)
            return false
        end
    end
    return true
end

function forwlossback(net, x, z)
    n = length(net)
    xx = Any[]
    for i=1:n
        x = forw(net[i], copy(x); seed=1)
        @assert isa(x, KUdense)
        push!(xx, x)
    end
    ll = (isa(net[n], Loss) ? loss(net[n],z) : 0)
    zz = Any[]
    for i=n:-1:1
        z = back(net[i], copy(z))
        @assert isa(z, KUdense)
        unshift!(zz, z)
    end
    return (ll, xx, zz)
end

function getloss(net, x, z)
    x = copy(x)
    isa(net[1],LogpLoss) && (x=forw(Logp(),x))
    isa(net[1],SoftLoss) && (x=KUdense(x.arr ./ sum(x.arr,1)))
    # @show x
    n = length(net)
    for i=1:n; x = forw(net[i], x; seed=1); @assert isa(x, KUdense); end
    loss(net[n], z)
end

function getparam(l1::Op)
    w1 = nothing
    for n in fieldnames(l1); isdefined(l1,n) && isa(l1.(n), KUparam) && (w1=l1.(n); break); end
    return w1
end

function gputest(cnet::MLP, x, z)
    global gnet0
    rval = true
    # Compare loss, y, dx, dw after forw and back:
    # info("gputest 1")
    # display(shownet(cnet));println("")
    (cl, cx, cz) = forwlossback(cnet, x, z)
    # info("gputest 2")
    # display(shownet(cnet));println("")
    gnet0 = gnet = Op[gpucopy(cnet)...]
    # info("gputest 3")
    # display(shownet(gnet));println("")
    # @show gnet
    # hnet = Op[cpucopy(gnet)...]
    # @show hnet
    # @assert isequal(cnet[1].w.arr, hnet[1].w.arr)
    # @show cnet
    xx,zz = gpucopy(x),gpucopy(z)
    # @show (xx,zz)
    (gl, gx, gz) = forwlossback(gnet, xx, zz)
    # info("gputest 4")
    # display(shownet(gnet));println("")
    isapprox(gl, cl; rtol=cbrt(eps(Float32)), atol=sqrt(eps(Float32))) || (warn("loss mismatch in $(map(typeof,cnet)): $gl != $cl"); rval=false)
    # info("gputest 5")
    for i=1:length(cnet)
        isapprox(cx[i], gx[i]) || (warn("y mismatch in $(typeof(cnet[i]))"); rval=false)
        isapprox(cz[i], gz[i]) || (warn("dx mismatch in $(typeof(cnet[i]))"); rval=false)
        cw = getparam(cnet[i]); gw = getparam(gnet[i])
        cw == nothing || isapprox(cw.diff, gw.diff) || (warn("dw mismatch in $(typeof(cnet[i]))"); rval=false)
    end

    # Compare w, dw after update:
    # info("gputest 6")
    setparam!(cnet; lr=0.1, l1reg=0.1, l2reg=0.1, adagrad=true, momentum=0.1, nesterov=0.1)
    # info("gputest 7")
    setparam!(gnet; lr=0.1, l1reg=0.1, l2reg=0.1, adagrad=true, momentum=0.1, nesterov=0.1)
    # info("gputest 8")
    update!(cnet)
    # info("gputest 9")
    update!(gnet)
    # info("gputest 10")
    for i=1:length(cnet)
        cw = getparam(cnet[i]); gw = getparam(gnet[i])
        cw == nothing && continue
        isapprox(cw.diff, gw.diff) || (warn("dw mismatch after update in $(typeof(cnet[i]))"); rval=false)
        isapprox(cw.arr, gw.arr) || (warn("w mismatch after update in $(typeof(cnet[i]))"); rval=false)
    end
    # info("gputest 11")

    return rval
end

function gradtest(net, x, z)
    rval = true
    (ll, yy, dx) = forwlossback(net, x, z)
    gradcheck(net, x, z, x, dx[1]) || (warn("gradtest failed for dx in $(typeof(net[1]))"); rval=false)
    w1 = getparam(net[1])
    w1 == nothing || gradcheck(net, x, z, w1.arr, w1.diff) || (warn("gradtest failed for dw in $(typeof(net[1]))"); rval=false)
    return rval
end

function iseq03(a,b)
    typeof(a)==typeof(b) || return false
    isa(a,Tuple) && return all(map(iseq03, a, b))
    isempty(fieldnames(a)) && return isequal(a,b)
    for n in fieldnames(a)
        in(n, (:x,:y,:dx,:dy,:xdrop)) && continue
        in(n, fieldnames(b)) || (warn("$n missing");return false)
        isdefined(a,n) || continue
        isdefined(b,n) || (warn("$n undefined");return false)
        iseq03(a.(n), b.(n)) || (warn("$n unequal"); return false)
    end
    for n in fieldnames(b)
        in(n, (:x,:y,:dx,:dy,:xdrop)) && continue
        in(n, fieldnames(a)) || (warn("$n missing");return false)
        isdefined(b,n) || continue
        isdefined(a,n) || (warn("$n undefined");return false)
    end
    return true
end

function filetest(net1)
    isa(net1[1], Pool) && (warn("Pooling layers cannot be saved to file yet"); return true)
    Knet.savenet("/tmp/knet.test", net1)
    net2 = Knet.loadnet("/tmp/knet.test")
    return all(map(iseq03, net1, net2))
end

function getnet{T<:Op}(F,S,L::Type{T})
    nd = length(S)
    nf = (nd==1 ? S[1] : div(prod(S),S[nd]))
    (nf>20) && in(L,(Logp,Soft,LogpLoss,SoftLoss,XentLoss)) && return nothing
    (nd!=4) && in(L,(Conv,Pool)) && return nothing
    C = (nd==1 ? S[1] : S[nd-1])
    l = ((L == Conv) ? Conv(rand(1:20), rand(1:min(S[1],S[2]))) :
         (L == Drop) ? Drop(rand()) :
         (L == Mmul) ? Mmul(rand(1:20)) :
         (L == Pool) ? Pool(rand(1:minimum(S))) : L())
    net = Op[]; push!(net, l)
    return (isa(l, Logp) ? push!(net, LogpLoss()) :
            isa(l, Soft) ? push!(net, SoftLoss()) :
            !isa(l, Loss) ? push!(net, QuadLoss()) : net)
end

function getx(F,S,L)
    x = KUdense(rand(F,S))
    return ((L == LogpLoss) ? forw(Logp(), x) :
            (L == SoftLoss) ? forw(Soft(), x) : x)
end

function getz(net, x)
    (net == nothing || x == nothing) && return nothing
    z = copy(rand!(forw(net, copy(x))))
    @assert isa(z, KUdense)
    L = typeof(net[end])
    return (in(L, (Logp,)) ? forw(Logp(), z) :
            in(L, (Soft, SoftLoss, LogpLoss, XentLoss)) ? forw(Soft(), z) : z)
end

function gettest(F,S,L)
    net = getnet(F,S,L)
    x = getx(F,S,L)
    z = getz(net, x)
    return (net, x, z)
end

function shownet(n::MLP)
    map(showlayer, n)
end

function showlayer(l::Op)
    ans = Any[]
    push!(ans,typeof(l))
    for n in fieldnames(l)
        isdefined(l,n) || continue
        f = l.(n)
        push!(ans,n)
        # push!(ans,typeof(f))
        if isa(f, KUdense)
            # @assert pointer(f) === pointer(f.arr) === pointer(f.ptr)
            push!(ans, pointer(f.arr))
        elseif isa(f, KUparam)
            push!(ans, pointer(f.arr))
        end
    end
    ans
end

function main(layers)
    global net0, x0, z0
    Knet.gpu(false)
    for F in (Float32,Float64)
        for D in 1:5
            S = tuple(rand(1:20,D)...)
            for L in layers
                (net, x, z) = gettest(F,S,L)
                net==nothing && continue  # combination not supported
                net0, x0, z0 = net, x, z
                @show (F, S, L)
                Knet.GPU && (@test gputest(net, x, z))
                gradtest(net, x, z)
                @test filetest(net)
            end
        end
    end
end

# Test each layer for: 1D-5D, gpu/cpu, float32/float64
layers = (Bias, Conv, Drop, Logp, LogpLoss, Mmul, PercLoss, Pool, QuadLoss, Relu, Sigm, Soft, SoftLoss, Tanh, XentLoss)
# These don't have cpu versions: Conv, Pool
# layers = (Bias, Drop, Logp, LogpLoss, Mmul, PercLoss, QuadLoss, Relu, Sigm, Soft, SoftLoss, Tanh, XentLoss)
# layers = (Conv,Pool)
main(layers)


text:
using Dates

include(joinpath("..", "src", "CairoPlot.jl"))
using .CairoPlot

# using Cairo
# include(joinpath("..", "src", "main.jl"))

x = today():Year(1):today()+Year(20)
y = rand(length(x))

xticks = today():Year(5):today()+Year(20)
yticks = 0.0:0.2:1.0

c = CairoPlot.crplot(x, y; xticks, yticks)
write_to_png(c,"cairoplot_dateexample.png")

nothing


text:
using NLopt, ForwardDiff

x0 = zeros(14)
gcfg = ForwardDiff.GradientConfig(nothing,x0)
g!(grad, x) = ForwardDiff.gradient!(grad, integral_function, x, gcfg)
lbs  = -5ones(14)
lbs[3:5] = zeros(3)
lbs[9] = 0.0
ubs = 20ones(14)
ubs[3:5] = ones(3)
ubs[9] = 1.0

function integral_function(x)
  a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032 = x
  b12 = (b22*(1 - c11 + b13*c11 - b13*c13))/(-1 + b23*c11 - b23*c13)
  b11 = 1 - b12 - b13
  b21 = - b22 - b23
  c12 = (-1 - b21*c11 - b23*c13)/b22
  a1 = 1 - a2 - a3
  c02 = A021
  c03 = A031 + A032
  G(A021,A031,A032,a1,a2,a3)
end

iter = 0
g_res = zeros(3)

function f(x::Vector, grad::Vector)
  if length(grad) > 0
    g!(grad,x)
  end
  integral = integral_function(x)
  global iter
  iter::Int += 1
  if mod(iter,1000) == 0
    global g_res
    g_res[1] = a2*B021 + a3*B031 + a3*B032 - 1
    g_res[2] = a2*B021^2 + a3*(B031 + B032)^2 - 3/2
    g_res[3] = A021*a2 + A031*a3 + A032*a3 - 1/2
    @show x
    @show g_res
    @show integral
  end
  integral
end

function G(A021,A031,A032,a1,a2,a3)
  tmps = zeros(Threads.nthreads())
  xs = -6:0.01:1
  ys = -4:0.01:4
  Threads.@threads for x in xs
    for y in ys
      z = x + im*y
      tmps[Threads.threadid()] += abs(1 + z*(a1 + a2 + A021*z*a2 + a3 + A032*z*a3 + (A031*z + A021*A032*z^2)*a3))<1
    end
  end
  tmp = sum(tmps)
  (tmp/(length(xs)*length(ys)))*(5*6)
end

function constraints!(res,x,grad)
  a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032 = x
  b12 = (b22*(1 - c11 + b13*c11 - b13*c13))/(-1 + b23*c11 - b23*c13)
  b11 = 1 - b12 - b13
  b21 = - b22 - b23
  c12 = (-1 - b21*c11 - b23*c13)/b22
  a1 = 1 - a2 - a3
  c02 = A021
  c03 = A031 + A032
  if length(grad)>0
    grad .= 0.0
    grad[1,1] = B021
    grad[2,1] = B031 + B032
    grad[12,1] = a2
    grad[13,1] = a3
    grad[14,1] = a3
    grad[1,2] = B021^2
    grad[2,2] = (B031 + B032)^2
    grad[12,2] = 2a2*B021
    grad[13,2] = 2a3*(B031+B032)
    grad[14,2] = 2a3*(B031+B032)
    grad[1,3] = A021
    grad[2,3] = A031 + A032
    grad[9,3] = a2
    grad[10,3] = a3
    grad[11,3] = a3
  end
  res[1] = a2*B021 + a3*B031 + a3*B032 - 1
  res[2] = a2*B021^2 + a3*(B031 + B032)^2 - 3/2
  res[3] = A021*a2 + A031*a3 + A032*a3 - 1/2
end

function c_constraints!(res,x,grad)
  a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032 = x
  b12 = (b22*(1 - c11 + b13*c11 - b13*c13))/(-1 + b23*c11 - b23*c13)
  b11 = 1 - b12 - b13
  b21 = - b22 - b23
  c12 = (-1 - b21*c11 - b23*c13)/b22
  a1 = 1 - a2 - a3
  c02 = A021
  c03 = A031 + A032
  res[1] = c03 - 1
  res[2] = c12 - 1
  res[3] = -c03
  res[4] = -c12
end

a2 = 1/6; a3 = 2/3; c01 = 0; c11 = 1; c13 = 0; b13 = 0;
b22 = 1; b23 = 0; A021 = 1; A031 = 1/4; A032 = 1/4;
B021 = 0; B031 = 1; B032 = 1/2
#x0 = [a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032] + rand(14)
x0 = rand(14)
x0[3:5] = max.(0,x0[3:5])
x0[3:5] = min.(1,x0[3:5])
x0[9] = max.(0,x0[9])
x0[9] = min.(1,x0[9])
opt = Opt(:LN_COBYLA, 14) #:LD_SLSQP, :LN_COBYLA (semi), :GN_ISRES, :LN_AUGLAG_EQ support equality constraints
lower_bounds!(opt, lbs)
upper_bounds!(opt, ubs)
max_objective!(opt, f)
equality_constraint!(opt::Opt,constraints!, [1e-12 for i in 1:3])
inequality_constraint!(opt::Opt,c_constraints!, [1e-12 for i in 1:4])
ftol_abs!(opt,1e-12)
xtol_rel!(opt,1e-12)
maxeval!(opt::Opt, Int(1e6))
maxf,maxx,ret = optimize(opt,x0)

#maxx_str = "x = [$(maxx[1]),$(maxx[2]),$(maxx[3]),$(maxx[4]),$(maxx[5]),$(maxx[6]),$(maxx[7]),$(maxx[8]),$(maxx[9]),$(maxx[10]),$(maxx[11]),$(maxx[12]),$(maxx[13]),$(maxx[14])]"
#println(maxx_str)

### Result

x = maxx
a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032 = x
b12 = (b22*(1 - c11 + b13*c11 - b13*c13))/(-1 + b23*c11 - b23*c13)
b11 = 1 - b12 - b13
b21 = - b22 - b23
c12 = (-1 - b21*c11 - b23*c13)/b22
a1 = 1 - a2 - a3
c02 = A021
c03 = A031 + A032
G(A021,A031,A032,a1,a2,a3)
g_res = zeros(3); g_res_grad = zeros(14,3); c_res=zeros(4)
constraints!(g_res, x, g_res_grad)
g_res
g_res_grad
c_constraints!(c_res, x, g_res_grad)
c_res

value_str = """
  a2 = $a2;
  a3 = $a3;
  c01 = $c01;
  c11 = $c11;
  c13 = $c13;
  b13 = $b13;
  b22 = $b22;
  b23 = $b23;
  A021 = $A021;
  A031 = $A031;
  A032 = $A032;
  B021 = $B021;
  B031 = $B031;
  B032 = $B032;
  b12 = $b12;
  b11 = $b11;
  b21 = $b21;
  c12 = $c12;
  a1 = $a1;
  c02 = $c02;
  c03 = $c03;
  """
println(value_str)
#=
##################
# Test Integral
G(1,13/56,-1/56,13/4,5/4,-7/2) # 12.063
G(1,1/4,1/4,1/6,1/6,2/3) # 9.0872

###################
# Test constraints

a2 = 1/6; a3 = 2/3; c01 = 0; c11 = 1; c13 = 0; b13 = 0;
b22 = 1; b23 = 0; A021 = 1; A031 = 1/4; A032 = 1/4;
B021 = 0; B031 = 1; B032 = 1/2
x = [a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032]

a2, a3, c01, c11, c13, b13, b22, b23, A021, A031, A032, B021, B031, B032 = x
b12 = (b22*(1 - c11 + b13*c11 - b13*c13))/(-1 + b23*c11 - b23*c13)
b11 = 1 - b12 - b13
b21 = - b22 - b23
c12 = (-1 - b21*c11 - b23*c13)/b22
a1 = 1 - a2 - a3
c02 = A021
c03 = A031 + A032

a1 ≈ 1/6
a2 ≈ 1/6
a3 == 2/3
c01 == 0
c02 == 1
c03 == 1/2
c11 == 1
c12 == 0
c13 == 0
A021 == 1
A031 == 1/4
A032 == 1/4
B021 == 0
B031 == 1
B032 == 1/2
b11 == 1
b12 == 0
b13 == 0
b21 == -1
b22 == 1
b23 == 0

g_res = zeros(3); g_res_grad = zeros(14,3); c_res=zeros(4)
constraints!(g_res, x, g_res_grad)
c_constraints!(c_res, x, g_res_grad)
g_res
g_res_grad
c_res
=#


text:
# This file is a part of JuliaFEM.
# License is MIT: see https://github.com/JuliaFEM/AbaqusReader.jl/blob/master/LICENSE

using AbaqusReader: abaqus_read_model
using AbaqusReader: parse_keyword

datadir = first(splitext(basename(@__FILE__)))

@testset "parse abaqus inp file to AbaqusModel" begin
    fn = joinpath(datadir, "cube_tet4.inp")
    model = abaqus_read_model(fn)

    @test length(model.properties) == 1
    section = first(model.properties)
    @test section.element_set == :CUBE
    @test section.material_name == :MAT

    @test haskey(model.materials, :MAT)
    material = model.materials[:MAT]
    @test isapprox(first(material.properties).E, 208.0e3)

    @test length(model.steps) == 1
    step = first(model.steps)
    @test length(step.boundary_conditions) == 2

    bc = step.boundary_conditions[1]
    @test bc.data[1] == [:SYM23, 1]
    @test bc.data[2] == [:SYM13, 2]

    load = step.boundary_conditions[2]
    @test load.data[1] == [:LOAD, :P, 1.00000]
end

@testset "parse keyword" begin
    k = parse_keyword("*SURFACE, NAME=TIE, SURFACE TO SURFACE")
    @test k.options[1] == ("NAME" => "TIE")
    @test k.options[2] == "SURFACE TO SURFACE"
    @test_throws(Exception, parse_keyword("*SURFACE, MOI=HEI=EI"))
end


text:
function instantiate_fplane_1(FT)
    coriolis = FPlane(FT, f=π)
    return coriolis.f == FT(π)
end

function instantiate_fplane_2(FT)
    coriolis = FPlane(FT, rotation_rate=2, latitude=30)
    return coriolis.f == FT(2)
end

function instantiate_betaplane_1(FT)
    coriolis = BetaPlane(FT, f₀=π, β=2π)
    return coriolis.f₀ == FT(π)
end

function instantiate_betaplane_2(FT)
    coriolis = BetaPlane(FT, latitude=70, radius=2π, rotation_rate=3π)
    return coriolis.f₀ == FT(6π * sind(70))
end

@testset "Coriolis" begin
    @info "Testing Coriolis..."

    @testset "Coriolis" begin
        for FT in float_types
            @test instantiate_fplane_1(FT)
            @test instantiate_fplane_2(FT)
            @test instantiate_betaplane_1(FT)
            @test instantiate_betaplane_2(FT)

            # Test that FPlane throws an ArgumentError
            @test_throws ArgumentError FPlane(FT, rotation_rate=7e-5)
            @test_throws ArgumentError FPlane(FT, f=1, latitude=40)
            @test_throws ArgumentError FPlane(FT, f=1, rotation_rate=7e-5, latitude=40)

            # Non-exhaustively test that BetaPlane throws an ArgumentError
            @test_throws ArgumentError BetaPlane(FT, f₀=1)
            @test_throws ArgumentError BetaPlane(FT, β=1)
            @test_throws ArgumentError BetaPlane(FT, f₀=1e-4, β=1e-11, latitude=70)
        end
    end
end


text:
using TensorCalculus
using Test

@testset "TensorCalculus" begin
    @testset "Tensor" begin
        include("tensor.jl")
    end

    @testset "Broadcasting" begin
        include("broadcast.jl")
    end

    @testset "Arithmetic" begin
        include("arithmetic.jl")
    end

    @testset "Algebra" begin
        include("algebra.jl")
    end
end


text:
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
#
#   Print function of the LaTeX backend.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #

# Low-level function to print the table using the LaTeX backend.
function _pt_latex(io, pinfo;
                   tf::LatexTableFormat = latex_default,
                   body_hlines::Vector{Int} = Int[],
                   cell_alignment::Dict{Tuple{Int,Int},Symbol} = Dict{Tuple{Int,Int},Symbol}(),
                   highlighters::Union{LatexHighlighter,Tuple} = (),
                   hlines::Union{Nothing,Symbol,AbstractVector} = nothing,
                   longtable_footer::Union{Nothing,AbstractString} = nothing,
                   noheader::Bool = false,
                   nosubheader::Bool = false,
                   row_number_alignment::Symbol = :r,
                   show_row_number::Bool = false,
                   table_type::Symbol = :tabular,
                   vlines::Union{Nothing,Symbol,AbstractVector} = nothing,
                   # Deprecated
                   formatter = nothing)

    @unpack_PrintInfo pinfo

    # We cannot use `@unpack_` here because it will overwrite `hlines` and
    # `vlines.`
    top_line       = tf.top_line
    header_line    = tf.header_line
    mid_line       = tf.mid_line
    bottom_line    = tf.bottom_line
    left_vline     = tf.left_vline
    mid_vline      = tf.mid_vline
    right_vline    = tf.right_vline
    header_envs    = tf.header_envs
    subheader_envs = tf.subheader_envs

    # Let's create a `IOBuffer` to write everything and then transfer to `io`.
    buf_io = IOBuffer()
    buf    = IOContext(buf_io)

    table_type ∉ [:tabular, :longtable] &&
    error("Unknown table type $table_type. Possible values are `:tabular` or `:longtable`.")

    table_env = table_type == :tabular ? "tabular" : "longtable"

    !noheader && num_cols != header_num_cols &&
    error("The header length must be equal to the number of columns.")

    # Additional processing necessary if the user wants to print the header.
    if !noheader
        # If the user do not want to print the sub-header but wants to print the
        # header, then just force the number of rows in header to be 1.
        if nosubheader
            # Now, `header` will be a view of the first line of the matrix that
            # has the header.
            header = @view header[1:header_num_rows:end]
            header_num_rows = 1
        end
    end

    # Make sure that `highlighters` is always a tuple.
    !(highlighters isa Tuple) && (highlighters = (highlighters,))

    # Get the string which is printed when `print` is called in each element of
    # the matrix. Notice that we must create only the matrix with the printed
    # rows and columns.
    header_str = Matrix{String}(undef, header_num_rows, num_printed_cols)
    data_str   = Matrix{AbstractString}(undef, num_printed_rows, num_printed_cols)

    @inbounds for i = 1:num_printed_cols
        # Index of the i-th printed column in `data`.
        ic = id_cols[i]

        if !noheader
            for j = 1:header_num_rows
                header_str[j,i] =
                    _str_latex_escaped(sprint(print, header[(ic-1)*header_num_rows + j];
                                              context = :compact => compact_printing))
            end
        end

        for j = 1:num_printed_rows
            # Index of the j-th printed row in `data`.
            jr = id_rows[j]

            # Apply the formatters.
            data_ij = isassigned(data,jr,ic) ? data[jr,ic] : undef

            for f in formatters
                data_ij = f(data_ij, jr, ic)
            end

            # Handle `nothing`, `missing`, and `undef`.
            if ismissing(data_ij)
                data_str_ij = "missing"
            elseif data_ij == nothing
                data_str_ij = "nothing"
            elseif data_ij == undef
                data_str_ij = "\\#undef"
            elseif data_ij isa Markdown.MD
                data_str_ij = replace(sprint(show, MIME("text/latex"), data_ij),"\n"=>"")
            else
                data_str_ij = sprint(print, data_ij;
                                     context = :compact => compact_printing)
            end

            # Check if the user wants to show only the first line.
            cell_first_line_only && (data_str_ij = split(data_str_ij, '\n')[1])

            data_str_ij_esc = _str_latex_escaped(data_str_ij)
            data_str[j,i]   = data_str_ij_esc
        end
    end

    # Compute where the horizontal and vertical lines must be drawn
    # --------------------------------------------------------------------------

    hlines == nothing && (hlines = tf.hlines)
    hlines = _process_hlines(hlines, body_hlines, num_printed_rows, noheader)

    # Process `vlines`.
    #
    # TODO: `num_printed_cols` must consider the row number.
    vlines == nothing && (vlines = tf.vlines)
    vlines = _process_vlines(vlines, num_printed_cols + show_row_number)

    # Variables to store information about indentation
    # ==========================================================================

    il = 0 # ......................................... Current indentation level
    ns = 2 # ........................ Number of spaces in each indentation level

    # Print LaTeX header
    # ==========================================================================

    if table_type == :tabular
        _aprintln(buf, "\\begin{table}", il, ns)
        il += 1
        length(title) > 0 && _aprintln(buf, "\\caption{$title}", il, ns)
    end

    _aprintln(buf,"""
              \\begin{$table_env}$(_latex_table_desc(id_cols,
                                                     alignment,
                                                     show_row_number,
                                                     row_number_alignment,
                                                     vlines,
                                                     left_vline,
                                                     mid_vline,
                                                     right_vline))""",
              il, ns)
    il += 1

    if table_type == :longtable
        length(title) > 0 && _aprintln(buf, "\\caption{$title}\\\\", il, ns)
    end

    0 ∈ hlines && _aprintln(buf, top_line, il, ns)

    # Data header
    # ==========================================================================

    # Header and sub-header texts
    # --------------------------------------------------------------------------

    if !noheader
        @inbounds @views for i = 1:header_num_rows
            # The text "Row" must appear only on the first line.
            if show_row_number
                _aprint(buf, il, ns)

                if i == 1
                    print(buf, _latex_envs("Row", header_envs))
                end

                print(buf, " & ")
            end

            for j = 1:num_printed_cols
                # Apply the alignment to this row.
                (!show_row_number && j == 1) && _aprint(buf, il, ns)

                # Index of the j-th printed column in `data`.
                jc = id_cols[j]

                # Configure the LaTeX environments for the header and
                # sub-headers.
                if i == 1
                    envs = header_envs
                else
                    envs = subheader_envs
                end

                header_str_ij = _latex_envs(header_str[i,j], envs)

                # Check the alignment of this cell.
                alignment_ij = header_alignment[jc]

                for f in header_cell_alignment
                    aux = f(header, i, jc)

                    if aux ∈ (:l, :c, :r, :L, :C, :R, :s, :S)
                        alignment_ij = aux
                        break
                    end
                end

                # If alignment is `:s`, then we must use the column alignment.
                alignment_ij ∈ (:s,:S) && (alignment_ij = alignment[jc])

                # Check if the alignment of the cell must be overridden.
                if alignment_ij != alignment[jc]
                    header_str_ij = _latex_apply_cell_alignment(header_str_ij,
                                                                alignment_ij, j,
                                                                num_printed_cols,
                                                                show_row_number,
                                                                vlines,
                                                                left_vline,
                                                                mid_vline,
                                                                right_vline)
                end

                print(buf, header_str_ij)

                j != num_printed_cols && print(buf, " & ")
            end

            print(buf, " \\\\")
            if (i == header_num_rows) && (1 ∈ hlines)
                print(buf, header_line)
            end
            println(buf, "")
        end
    end

    # If we are using `longtable`, then we must mark the end of header and also
    # create the footer.
    if table_type == :longtable
        _aprintln(buf, "\\endhead", il, ns)
        _aprintln(buf, bottom_line, il, ns)

        # Check if the user wants a text on the footer.
        if longtable_footer != nothing
            lvline =            0 ∈ vlines ? left_vline : ""
            rvline = id_cols[end] ∈ vlines ? right_vline : ""

            env = "multicolumn{" * string(num_printed_cols) * "}" * "{r}"

            _aprintln(buf, _latex_envs(longtable_footer, env) * "\\\\", il, ns)
            _aprintln(buf, bottom_line, il, ns)
        end

        _aprintln(buf, "\\endfoot", il, ns)
        _aprintln(buf, "\\endlastfoot", il, ns)
    end

    # Data
    # ==========================================================================

    @inbounds @views for i = 1:num_printed_rows
        ir = id_rows[i]

        if show_row_number
            _aprint(buf, string(ir) * " & ", il, ns)
        end

        for j = 1:num_printed_cols
            # Apply the alignment to this row.
            (!show_row_number && j == 1) && _aprint(buf, il, ns)

            jc = id_cols[j]

            # If we have highlighters defined, then we need to verify if this
            # data should be highlight.
            data_str_ij = data_str[i,j]

            for h in highlighters
                if h.f(_getdata(data), ir, jc)
                    data_str_ij = h.fd(_getdata(data), i, j, data_str[i,j])
                    break
                end
            end

            # Check the alignment of this cell.
            alignment_ij = alignment[jc]

            for f in cell_alignment
                aux = f(_getdata(data), ir, jc)

                if aux ∈ [:l, :c, :r, :L, :C, :R]
                    alignment_ij = aux
                    break
                end
            end

            # Check if the alignment of the cell must be overridden.
            if alignment_ij != alignment[jc]
                data_str_ij = _latex_apply_cell_alignment(data_str_ij,
                                                          alignment_ij, j,
                                                          num_printed_cols,
                                                          show_row_number,
                                                          vlines, left_vline,
                                                          mid_vline,
                                                          right_vline)
            end

            print(buf, data_str_ij)
            j != num_printed_cols && print(buf, " & ")
        end

        print(buf, " \\\\")

        # Check if the user wants a horizontal line here.
        if (i+!noheader) ∈ hlines
            if i != num_printed_rows
                # Check if we must draw a horizontal line here.
                print(buf, mid_line)
            else
                print(buf, bottom_line)
            end
        end
        println(buf, "")
    end

    # Print LaTeX footer
    # ==========================================================================

    il -= 1
    _aprintln(buf, "\\end{$table_env}", il, ns)

    if table_type == :tabular
        il -= 1
        _aprintln(buf, "\\end{table}", il, ns)
    end

    # Print the buffer into the io.
    # ==========================================================================

    print(io, String(take!(buf_io)))

    return nothing
end


text:
#   This file is part of Maxima.jl. It is licensed under the MIT license
#   Copyright (c) 2016 Nathan Smith
import REPL

Reset(;args::Cmd=``) = (kill(ms); Load(args=args))
__init__() = (Load(); atexit(() -> kill(ms)))

function Load(;args::Cmd=``)
    try
        Sys.isunix() ? read(`maxima --version`, String) : read(`maxima.bat --version`, String)
    catch err
        error("Looks like Maxima is either not installed or not in the path")
    end

    # Server setup

    global ms = MaximaSession(args=args)	# Spin up a Maxima session

    # REPL setup
    repl_active = isdefined(Base, :active_repl)	# Is an active repl defined?
    interactive = isinteractive()				# In interactive mode?

    if repl_active && interactive
        repl_init(Base.active_repl)
    else # package is loaded from ~/.juliarc.jl
        atreplinit() do repl # check if OhMyREPL is loaded
            !isdefined(Main, :OhMyREPL) && (repl.interface = REPL.setup_interface(repl))
            repl_init(Base.active_repl)
        end
    end

    mcall("1") # clears out any initial data
    return nothing
end


text:
# SPDX-License-Identifier: X11
# 2020-09-13 (Live)
# Not (100pt)

print(parse(Int, readline()) == 1 ? 0 : 1)

max_stars_repo_name: jthvai/competitive-programming

max_issues_repo_name: jthvai/competitive-programming

max_forks_repo_name: jthvai/competitive-programming


text:
# Contiguous on-device arrays

export CuDeviceArray, CuDeviceVector, CuDeviceMatrix, ldg

## construction

"""
    CuDeviceArray(dims, ptr)
    CuDeviceArray{T}(dims, ptr)
    CuDeviceArray{T,N}(dims, ptr)
    CuDeviceArray{T,N,A}(dims, ptr)

Construct an `N`-dimensional dense CUDA device array with element type `T` wrapping a
pointer, where `N` is determined from the length of `dims` and `T` is determined from the
type of `ptr`. `dims` may be a single scalar, or a tuple of integers corresponding to the
lengths in each dimension). If the rank `N` is supplied explicitly as in `Array{T,N}(dims)`,
then it must match the length of `dims`. The same applies to the element type `T`, which
should match the type of the pointer `ptr`.
"""
CuDeviceArray

# NOTE: we can't support the typical `tuple or series of integer` style construction,
#       because we're currently requiring a trailing pointer argument.

struct CuDeviceArray{T,N,A} <: AbstractArray{T,N}
    shape::Dims{N}
    ptr::LLVMPtr{T,A}

    # inner constructors, fully parameterized, exact types (ie. Int not <:Integer)
    CuDeviceArray{T,N,A}(shape::Dims{N}, ptr::LLVMPtr{T,A}) where {T,A,N} = new(shape,ptr)
end

const CuDeviceVector = CuDeviceArray{T,1,A} where {T,A}
const CuDeviceMatrix = CuDeviceArray{T,2,A} where {T,A}

# outer constructors, non-parameterized
CuDeviceArray(dims::NTuple{N,<:Integer}, p::LLVMPtr{T,A})                where {T,A,N} = CuDeviceArray{T,N,A}(dims, p)
CuDeviceArray(len::Integer,              p::LLVMPtr{T,A})                where {T,A}   = CuDeviceVector{T,A}((len,), p)

# outer constructors, partially parameterized
CuDeviceArray{T}(dims::NTuple{N,<:Integer},   p::LLVMPtr{T,A}) where {T,A,N} = CuDeviceArray{T,N,A}(dims, p)
CuDeviceArray{T}(len::Integer,                p::LLVMPtr{T,A}) where {T,A}   = CuDeviceVector{T,A}((len,), p)
CuDeviceArray{T,N}(dims::NTuple{N,<:Integer}, p::LLVMPtr{T,A}) where {T,A,N} = CuDeviceArray{T,N,A}(dims, p)
CuDeviceVector{T}(len::Integer,               p::LLVMPtr{T,A}) where {T,A}   = CuDeviceVector{T,A}((len,), p)

# outer constructors, fully parameterized
CuDeviceArray{T,N,A}(dims::NTuple{N,<:Integer}, p::LLVMPtr{T,A}) where {T,A,N} = CuDeviceArray{T,N,A}(Int.(dims), p)
CuDeviceVector{T,A}(len::Integer,               p::LLVMPtr{T,A}) where {T,A}   = CuDeviceVector{T,A}((Int(len),), p)

## getters

Base.pointer(a::CuDeviceArray) = a.ptr
Base.pointer(a::CuDeviceArray, i::Integer) =
    pointer(a) + (i - 1) * Base.elsize(a)

Base.elsize(::Type{<:CuDeviceArray{T}}) where {T} = sizeof(T)
Base.size(g::CuDeviceArray) = g.shape
Base.length(g::CuDeviceArray) = prod(g.shape)

Base.sizeof(x::CuDeviceArray) = Base.elsize(x) * length(x)

## conversions

Base.unsafe_convert(::Type{LLVMPtr{T,A}}, a::CuDeviceArray{T,N,A}) where {T,A,N} = pointer(a)

## indexing intrinsics

# NOTE: these intrinsics are now implemented using plain and simple pointer operations;
#       when adding support for isbits union arrays we will need to implement that here.

# TODO: arrays as allocated by the CUDA APIs are 256-byte aligned. we should keep track of
#       this information, because it enables optimizations like Load Store Vectorization
#       (cfr. shared memory and its wider-than-datatype alignment)

@inline function arrayref(A::CuDeviceArray{T}, index::Int) where {T}
    @boundscheck checkbounds(A, index)
    align = Base.datatype_alignment(T)
    unsafe_load(pointer(A), index, Val(align))
end

@inline function arrayset(A::CuDeviceArray{T}, x::T, index::Int) where {T}
    @boundscheck checkbounds(A, index)
    align = Base.datatype_alignment(T)
    unsafe_store!(pointer(A), x, index, Val(align))
    return A
end

@inline function const_arrayref(A::CuDeviceArray{T}, index::Int) where {T}
    @boundscheck checkbounds(A, index)
    align = Base.datatype_alignment(T)
    unsafe_cached_load(pointer(A), index, Val(align))
end

## indexing

Base.@propagate_inbounds Base.getindex(A::CuDeviceArray{T}, i1::Int) where {T} =
    arrayref(A, i1)
Base.@propagate_inbounds Base.setindex!(A::CuDeviceArray{T}, x, i1::Int) where {T} =
    arrayset(A, convert(T,x)::T, i1)

Base.IndexStyle(::Type{<:CuDeviceArray}) = Base.IndexLinear()

## const indexing

"""
    Const(A::CuDeviceArray)

Mark a CuDeviceArray as constant/read-only. The invariant guaranteed is that you will not
modify an CuDeviceArray for the duration of the current kernel.

This API can only be used on devices with compute capability 3.5 or higher.

!!! warning
    Experimental API. Subject to change without deprecation.
"""
struct Const{T,N,AS} <: DenseArray{T,N}
    a::CuDeviceArray{T,N,AS}
end
Base.Experimental.Const(A::CuDeviceArray) = Const(A)

Base.IndexStyle(::Type{<:Const}) = IndexLinear()
Base.size(C::Const) = size(C.a)
Base.axes(C::Const) = axes(C.a)
Base.@propagate_inbounds Base.getindex(A::Const, i1::Int) = const_arrayref(A.a, i1)

# deprecated
Base.@propagate_inbounds ldg(A::CuDeviceArray, i1::Integer) = const_arrayref(A, Int(i1))

## other

Base.show(io::IO, a::CuDeviceVector) =
    print(io, "$(length(a))-element device array at $(pointer(a))")
Base.show(io::IO, a::CuDeviceArray) =
    print(io, "$(join(a.shape, '×')) device array at $(pointer(a))")

Base.show(io::IO, mime::MIME"text/plain", a::CuDeviceArray) = show(io, a)

@inline function Base.iterate(A::CuDeviceArray, i=1)
    if (i % UInt) - 1 < length(A)
        (@inbounds A[i], i + 1)
    else
        nothing
    end
end


text:
using LightGraphs, Plots, StatsPlots, Distributions, GraphRecipes, GraphIO, ParserCombinator, JSON
using GraphIO: GML
gr()
theme(:juno)

N = 50
k = 1

function gen_graph(N, k)
    G = barabasi_albert(N, k)

    path = a_star(G, 1, N)

    uni = Set{Int}()
    for e ∈ path
        push!(uni, e.src)
        push!(uni, e.dst)
    end

    labels = Dict(i => i in uni ? 1 : 0 for i in 1:N)
    return G, labels
end

function plot_syn_graph(g, labels)
    graphplot(g, nodecolor=map((x) -> labels[x] == 0 ? :green : :blue, 1:length(labels)))
end

function write_to_file(g, labels, n)
    open("/gpfs_home/spate116/singhlab/GCN_Integration/scripts/BI/path/data/syn_graph_labels_$(n).json", "w") do io
        write(io, JSON.json(labels))
    end;

    open("/gpfs_home/spate116/singhlab/GCN_Integration/scripts/BI/path/data/syn_graph_$(n).gml", "w") do io
        GML.savegml(io, g)
    end;
end

G, labels = gen_graph(N, k)
write_to_file(G, labels, 0)

for i in 1:length(labels)
    if labels[i] == 1
        print(i)
        print('\n')
    end
end


text:
using Test
using MieSeries

@testset "Test Hankel" begin
    x = 4
    nu = 0:4
    sh, sh_d, sh_dd = MieSeries.sphericalhankel2_and_derivatives(nu, x)
    @test sh[1] ≈ -0.18920062382698205-0.16341090521590299im  # order 0
    @test sh[2] ≈ 0.11611074925915747-0.2300533501309578im
    @test sh[3] ≈ 0.27628368577135015-0.009129107382315343im
    @test sh[4] ≈ 0.22924385795503022+0.21864196590306362im
    @test sh[5] ≈ 0.12489306564995283+0.39175254771267665im
    @test sh_d[1] ≈ -0.11611074925915747+0.2300533501309578im  # order 0
    @test sh_d[2] ≈ -0.2472559984565608-0.04838423015042409im
    @test sh_d[3] ≈ -0.09110201506935513-0.22320651959422128im
    @test sh_d[4] ≈ 0.047039827816319935-0.22777107328537896im
    @test sh_d[5] ≈ 0.0731275258925892-0.2710487187377822im
    @test sh_dd[1] ≈ 0.247255998456560+0.0483842301504240im  # order 0
    @test sh_dd[2] ≈ 0.022031093626517+0.2254887964398001im
    @test sh_dd[3] ≈ -0.127126296072416+0.1173089519110577im
    @test sh_dd[4] ≈ -0.080830878396917+0.0592250451669235im
    @test sh_dd[5] ≈ -0.005340496533806+0.2334624962970602im

    for (index,i) in enumerate(nu)
        @test MieSeries.sphericalhankel2_and_derivatives(i, x) == (sh[index], sh_d[index], sh_dd[index])
    end
    @test MieSeries.sphericalhankel1_and_derivatives(nu, x) == conj.((sh,sh_d,sh_dd))

    sh2, sh_d2, sh_dd2 = MieSeries.sphericalhankel2_and_derivatives(1:4, x)
    @test sh2 == sh[2:end]
    @test sh_d2 == sh_d[2:end]
    @test sh_dd2 == sh_dd[2:end]

    # Spherical Bessel differential equation
    nu = 0:50
    sh, sh_d, sh_dd = MieSeries.sphericalhankel2_and_derivatives(nu, x)
    for n in nu
        index = n+1
        @test sh_dd[index] ≈ -(2*x*sh_d[index] + (x^2 - n*(n+1))*sh[index])/x^2
    end
end

@testset "Test Riccati-Bessel Hankel" begin
    x = 3.3
    nu = 0:2
    rh, rh_d, rh_dd = MieSeries.riccatihankel2_and_derivatives(nu, x)
    @test rh[1] ≈ -0.15774569-0.98747977im  # order 0
    @test rh[2] ≈ 0.93967804-0.45698199im
    @test rh[3] ≈ 1.01199846+0.5720416im
    @test rh_d[1] ≈ -0.98747977+0.15774569im  # order 0
    @test rh_d[2] ≈ -0.44249662-0.84900038im
    @test rh_d[3] ≈ 0.32634564-0.80367387im

    # Riccati-Bessel differential equation
    nu = 0:50
    rh, rh_d, rh_dd = MieSeries.riccatihankel2_and_derivatives(nu, x)
    for n in nu
        index = n+1
        @test rh_dd[index] ≈ -(1-n*(n+1)/x^2)*rh[index]
    end
end

max_stars_repo_name: Riarrieta/MieSeries

max_issues_repo_name: Riarrieta/MieSeries

max_forks_repo_name: Riarrieta/MieSeries


text:
"""
`FiniteDifferenceBVPProblem` - Initalize the data structure
"""
function FiniteDifferenceBVPProblem(a, b, n, f)
    return FiniteDifferenceBVPProblem(a, b, n, f, DirichletBC(0), DirichletBC(0))
end

function FiniteDifferenceBVPProblem(a, b, n, f, bcl::TBCL, bcr::TBCR) where {TBCL<:DirichletBC, TBCR<:DirichletBC}
    Δx = (b-a)/(n+1);
    x = LinRange(a, b, n+2)[2:end-1];

    # allocate the matrix
    A = zeros(n, n);
    rhs = zeros(n);

    return FiniteDifferenceBVPProblem(a, b, Δx, n, x, f, A, rhs,
        bcl, bcr)
end

function FiniteDifferenceBVPProblem(a, b, n, f, bcl::TBCL, bcr::TBCR) where {TBCL<:NeumannBC, TBCR<:NeumannBC}
    Δx = (b-a)/(n+1);
    x = LinRange(a, b, n+2);

    # allocate the matrix
    A = zeros(n+2, n+2);
    rhs = zeros(n+2);

    return FiniteDifferenceBVPProblem(a, b, Δx, n, x, f, A, rhs,
        bcl, bcr)
end

function FiniteDifferenceBVPProblem(a, b, n, f, bcl::TBCL, bcr::TBCR) where {TBCL<:PeriodicBC, TBCR<:PeriodicBC}
    Δx = (b-a)/(n+1);
    x = LinRange(a, b, n+2)[1:end-1];

    # allocate the matrix
    A = zeros(n+1, n+1);
    rhs = zeros(n+1);

    return FiniteDifferenceBVPProblem(a, b, Δx, n, x, f, A, rhs,
        bcl, bcr)
end


text:
using PSPipePlanckRender
using Test

@testset "PSPipePlanckRender.jl" begin
    # Write your tests here.
end

max_stars_repo_name: xzackli/PSPipePlanckRender

max_issues_repo_name: xzackli/PSPipePlanckRender

max_forks_repo_name: xzackli/PSPipePlanckRender


text:
#try with sparse
const STEP = 1
const SPIKE = 2
const SLOPE = 3
const SIN = 4
const COS = 5
const TOTALCOMPONENTS = 5

type iterator
	obs::Int

	components::Vector{Int}
    #??
	elements::Vector{Vector{Int}}#Vector{UnitRange{Int}}#Union{Vector{UnitRange{Int}},Vector{Vector{Int}}}

	#iterating_elements::Union{Vector{UnitRange{Int}},Vector{Vector{Int}}}
	# inicilaiza igual ao de cima
	# iteraness em vez do de cima
	# depois do primeiro lasso modidifica esse

	nelements::Vector{Int}

	isrange::Vector{Bool}

	ttelements::Int

	maxIter::Int
end

type dataCD

	#y0::Float64

	fs::Vector{Float64}
	fc::Vector{Float64}

	μ::Vector{Vector{Float64}}
	σ::Vector{Vector{Float64}}

	function dataCD()
		temp1 = Vector{Float64}[]
		temp2 = Vector{Float64}[]
		for i in 1:TOTALCOMPONENTS
			push!(temp1,zeros(0))
			push!(temp2,zeros(0))
		end
		new(zeros(0),
			zeros(0),
			temp1,
			temp2)
	end
end

max_stars_repo_name: joaquimg/ComponentFilterCD

max_issues_repo_name: joaquimg/ComponentFilterCD

max_forks_repo_name: joaquimg/ComponentFilterCD


text:
mutable struct Neighbourhood
    x::Matrix{Float64}
    f::Float64
    obj::Vector{Float64}
    infeas::Vector{Float64}
    ol::Vector{Float64}
    iu::Float64
    Neighbourhood(x, f, obj, infeas, ol, iu) = new(x, f, obj, infeas, ol, iu)
    Neighbourhood() = new(
        Matrix{Float64}(undef, 0, 0),
        Inf,
        Float64[],
        Float64[],
        Float64[],
        zero(Float64),
    )
end


text:
# MIT License
#
# Copyright (c) 2020 Joshua E Gentry

# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in all
# copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
# SOFTWARE.
using Test

# =========================================================================================
# tess_pipeline_hocr() => Filename
@testset "tess_pipeline_hocr(fontInfo = $font) => Filename" for font in [ true, false ]
    local inst       = TessInst("eng", datadir)
    local pipeline   = TessPipeline(inst)
    local output, io = mktemp(;cleanup=false); close(io)

    @test tess_pipeline_hocr(pipeline, output, font) == true

    @test (tess_run_pipeline(pipeline) do add
        @test add(pix_with([Line("This is image one. ®")]), 72) == true
        @test add(pix_with([Line("This is image two.")]), 72) == true
        @test add(pix_with([Line("This is image three.")]), 72) == true
        @test add(pix_with([Line("This is image four.")]), 72) == true
    end) == true

    @test isfile(output) == true
    @test filesize(output) > 0

    local text  = read(output, String)
    local lines = split(text, "\n"; keepempty=false)
    @test length(lines) == 62

    @test strip(lines[1]) == "<?xml version=\"1.0\" encoding=\"UTF-8\"?>"
    if font == true
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf ocrp_lang ocrp_dir ocrp_font ocrp_fsize'/>"
    else
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf'/>"
    end
    @test strip(lines[20]) == "®"
    @test strip(lines[31]) == "image"
    @test strip(lines[41]) == "This"
    @test strip(lines[51]) == """
"""
    @test strip(lines[61]) == ""

    rm(output)
end

# =========================================================================================
# tess_pipeline_hocr() => String
@testset "tess_pipeline_hocr(fontInfo = $font) => String" for font in [ true, false ]
    local inst     = TessInst("eng", datadir)
    local pipeline = TessPipeline(inst)
    local output   = tess_pipeline_hocr(pipeline, font)

    @test is_available(output) == false

    @test (tess_run_pipeline(pipeline) do add
        @test add(pix_with([Line("This is image one. ®")]), 72) == true
        @test add(pix_with([Line("This is image two.")]), 72) == true
        @test add(pix_with([Line("This is image three.")]), 72) == true
        @test add(pix_with([Line("This is image four.")]), 72) == true
    end) == true

    @test is_available(output) == true
    @test output[] != nothing

    local lines = split(output[], "\n"; keepempty=false)
    @test length(lines) == 62

    @test strip(lines[1]) == "<?xml version=\"1.0\" encoding=\"UTF-8\"?>"
    if font == true
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf ocrp_lang ocrp_dir ocrp_font ocrp_fsize'/>"
    else
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf'/>"
    end
    @test strip(lines[20]) == "®"
    @test strip(lines[31]) == "image"
    @test strip(lines[41]) == "This"
    @test strip(lines[51]) == """
"""
    @test strip(lines[61]) == ""
end

# =========================================================================================
# tess_pipeline_hocr() => Callback
@testset "tess_pipeline_hocr(fontInfo = $font) => Callback" for font in [ true, false ]
    local inst     = TessInst("eng", datadir)
    local pipeline = TessPipeline(inst)
    local lines    = Vector{String}()

    tess_pipeline_hocr(pipeline, font) do line
        push!(lines, line)
    end

    @test (tess_run_pipeline(pipeline) do add
        @test add(pix_with([Line("This is image one. ®")]), 72) == true
        @test add(pix_with([Line("This is image two.")]), 72) == true
        @test add(pix_with([Line("This is image three.")]), 72) == true
        @test add(pix_with([Line("This is image four.")]), 72) == true
    end) == true

    @test length(lines) == 62

    @test strip(lines[1]) == "<?xml version=\"1.0\" encoding=\"UTF-8\"?>"
    if font == true
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf ocrp_lang ocrp_dir ocrp_font ocrp_fsize'/>"
    else
        @test strip(lines[9]) == "<meta name='ocr-capabilities' content='ocr_page ocr_carea ocr_par ocr_line ocrx_word ocrp_wconf'/>"
    end
    @test strip(lines[20]) == "®"
    @test strip(lines[31]) == "image"
    @test strip(lines[41]) == "This"
    @test strip(lines[51]) == """
"""
    @test strip(lines[61]) == ""
end


text:
# ps3000_get_streaming_values_no_aggregation.jl

export ps3000_get_streaming_values_no_aggregation!, StreamingBufferNoAggregation

immutable StreamingBufferNoAggregation
  a :: Array{Int16,1}
  b :: Array{Int16,1}
  c :: Array{Int16,1}
  d :: Array{Int16,1}
  StreamingBufferNoAggregation(length) = new(Array(Int16,length),Array(Int16,length),
                                Array(Int16,length),Array(Int16,length))
end

function ps3000_get_streaming_values_no_aggregation!(handle        :: Int16,
                                                      st           :: Float64,
                                                      buffer       :: StreamingBufferNoAggregation,
                                                      no_of_values)
  start_time = Ref{Float64}(st)   # start, relative to trigger, time in ns
  overflow = Ref{Int16}           # bitfield, see doc
  triggered_at = Ref{UInt32}      # index trigger occured at
  triggered = Ref{UInt32}         # was there a valid trigger?
    number_written_to_buffer = ccall((:ps3000_get_streaming_values_no_aggregation, ps3000driver),UInt32,
    (Int16, Ref{Cdouble},
    Ptr{Int16},
    Ptr{Int16},
    Ptr{Int16},
    Ptr{Int16},
    Ref{Int16}, Ref{UInt32}, Ref{Int16}, UInt32),
    handle, start_time,
    buffer.a,
    buffer.b,
    buffer.c,
    buffer.d,
    overflow, triggered_at, triggered, no_of_values)
  if number_written_to_buffer == 0
    error("ps3000 Error: Parameter Out Of Range")
  end
  return (number_written_to_buffer,overflow[],triggered_at[],triggered[]==0?false:true)
end


text:
function interpolate_state(mdp::PedCarMDP, s::PedCarMDPState)
    # interpolate s in the grid
    vspace = get_car_vspace(mdp.env, mdp.vel_res)
    v_ped_space = get_ped_vspace(mdp.env, mdp.vel_ped_res)
    itp_ped, itp_ped_w = interpolate_pedestrian(get_ped_mdp(mdp), s.ped)
    itp_car, itp_car_w = interpolate_state(get_car_mdp(mdp, s.route), s.car)
    itp_ego, itp_ego_w = interpolate_state(mdp, s.ego)
    itp_states = Vector{PedCarMDPState}(length(itp_ego)*length(itp_car)*length(itp_ped))
    itp_w = Vector{Float64}(length(itp_states))
    l = 1
    for (i, p) in enumerate(itp_ped)
        for (j, c) in enumerate(itp_car)
            for (k, e) in enumerate(itp_ego)
                itp_states[l] = PedCarMDPState(crash(mdp, s), e, p, c, s.route)
                itp_w[l] = itp_ped_w[i]*itp_car_w[j]*itp_ego_w[k]
                l += 1
            end
        end
    end
    @assert sum(itp_w) ≈ 1.
    return itp_states, itp_w
end

function get_mdp_state(mdp::PedCarMDP, pomdp::UrbanPOMDP, s::Scene, ped_id, car_id)
    car_i = findfirst(car_id, s)
    car = Vehicle(get_off_the_grid(mdp), mdp.car_type, car_id)
    if car_i != nothing
        car = s[car_i]
    end
    ped_i = findfirst(ped_id, s)
    ped = Vehicle(get_off_the_grid(mdp), mdp.ped_type, ped_id)
    if ped_i != nothing
        ped = s[ped_i]
    end
    ego = get_ego(s)
    # find route
    sroute = nothing
    if haskey(pomdp.models, car_id) && car_i != nothing
        # find the exact route from the list of routes
        curr_route = [l.tag for l in pomdp.models[car_id].navigator.route]
        for route in get_car_routes(mdp.env)
            tags = intersect(Set(curr_route), Set(route))
            if length(tags) >= 2
                sroute = SVector{2, LaneTag}(route[1], route[end])
            elseif length(curr_route) == 1 && curr_route[1] ∈ route
                sroute = SVector{2, LaneTag}(route[1], route[end])
            end
        end
        if sroute == nothing
            println(curr_route)
        end
    else
        sroute = SVector{2, LaneTag}(LaneTag(0,0), LaneTag(0,0))
    end
    e_state = VehicleState(ego.state.posG, car_roadway(mdp.env), ego.state.v)
    p_state = VehicleState(ped.state.posG, mdp.env.ped_roadway, ped.state.v)
    c_state = VehicleState(car.state.posG, car_roadway(mdp.env), car.state.v)
    return PedCarMDPState(is_colliding(ego, car), e_state, p_state, c_state, sroute)
end

function get_car_mdp(mdp::PedCarMDP, route::SVector{2, LaneTag})
    # return a corresponding CarMDP
    if route == SVector(LaneTag(0,0), LaneTag(0,0))
       model = mdp.car_models[SVector(LaneTag(1,1), LaneTag(2,1))]
    else
       model = mdp.car_models[route]
    end
    return CarMDP(mdp.env, mdp.ego_type, mdp.car_type, model, mdp.max_acc,
                  mdp.pos_res, mdp.vel_res, mdp.ego_start, mdp.ego_goal, mdp.off_grid,
                  mdp.ΔT, mdp.car_birth, mdp.collision_cost, mdp.action_cost,
                  mdp.goal_reward, mdp.γ, mdp._ped_grid, mdp._car_grid,
                  mdp._l_grid, mdp._v_grid, mdp._collision_checker)
end

function get_ped_mdp(mdp::PedCarMDP)
    # return a corresponding CarMDP
    return PedMDP(mdp.env, mdp.ego_type, mdp.ped_type, mdp.ped_model, mdp.max_acc,
                  mdp.pos_res, mdp.vel_res, mdp.vel_ped_res, mdp.ego_start, mdp.ego_goal, mdp.off_grid,
                  mdp.ΔT, mdp.ped_birth, mdp.collision_cost, mdp.action_cost,
                  mdp.goal_reward, mdp.γ, mdp._ped_grid, mdp._car_grid,
                  mdp._l_grid, mdp._v_grid, mdp._collision_checker)
end

max_stars_repo_name: MaximeBouton/AutomotivePOMDPs

max_issues_repo_name: MaximeBouton/AutomotivePOMDPs

max_forks_repo_name: MaximeBouton/AutomotivePOMDPs


text:
"""
Transform the integral problem based on the supplied properties.

For example, if the integrand is known to be singular, then the domain may be
split into several components such that there are no singularities in the
interior.
"""
process_properties(qs, integrand, domain, measure, properties...) =
    (integrand, domain, measure, properties...)

#################
## Singularities
#################

filter_singularities() = ()
filter_singularities(property::Singularity) = ()
filter_singularities(property::Property) = (property,)
filter_singularities(property::Singularity, properties...) =
    filter_singularities(properties...)
filter_singularities(property::Property, properties...) =
    (property, filter_singularities(properties...)...)

property_splits_domain(property) = false
property_splits_domain(property::PointSingularity) = true
property_splits_domain(property::SingularDiagonal) = true

property_splits_domain() = false
property_splits_domain(property1, property2, properties...) =
    property_splits_domain(property1) || property_splits_domain(property2, properties...)

property_split_domain(domain) = (domain,)
property_split_domain(domain, property) = splitdomain_sing(property, domain)

function property_split_domain(domain, property, properties...)
    domains = Domain[]
    push!(domains, domain)
    property_split_multiple(domains, property, properties...)
end

property_split_multiple(domains) = domains
function property_split_multiple(domains_to_split, property, properties...)
    if property_splits_domain(property)
        domains = Domain[]
        for domain_to_split in domains_to_split
            for d in property_split_domain(domain_to_split, property)
                push!(domains, d)
            end
        end
        property_split_multiple(domains, properties...)
    else
        property_split_multiple(domains_to_split, properties...)
    end
end

"""
Split the domain into a list of domains, such that any singularity lies
on the boundary of one or more of the constituting parts.

Care is taken to avoid a singularity in the interior due to roundoff errors.
"""
splitdomain_sing(sing, domain) = (domain,)

splitdomain_sing(sing::PointSingularity, domain) = splitdomain_point(point(sing), domain)

function splitdomain_point(point, domain::MappedDomain)
    m = inverse_map(domain)
    map_domain.(Ref(m), splitdomain_point(m(point), superdomain(domain)))
end

function splitdomain_point(point, domain::Domain)
    if point ∈ domain
        @warn("Don't know how to split domain $(domain) at point $(point)")
        (domain,)
    else
        (domain,)
    end
end

"Split the domain according to a singularity at the point `x`."
function splitdomain_point(x, domain::AbstractInterval)
    T = promote_type(typeof(x),prectype(domain))
    tol = tolerance(T)
    a, b = extrema(domain)
    if (x > a+tol) && (x < b-tol)
        # x is in [a,b] and sufficiently far away from the boundary
        (a..x-tol, x+tol..b)
    elseif (x < a-tol) || (x > b+tol)
        # x is not in [a,b] and sufficiently far away from the boundary
        (domain,)
    else
        # x is close to a or b
        if abs(x-a) <= tol
            (max(a,x)+tol..b,)
        else
            (a..min(b,x)-tol,)
        end
    end
end

splitdomain_point(x, domain::ProductDomain) = splitdomain_point(x, domain, components(domain)...)

function splitdomain_point(x, domain::ProductDomain, domain1, domain2)
    domains1 = splitdomain_point(x[1], domain1)
    domains2 = splitdomain_point(x[2], domain2)
    [d1 × d2 for d1 in domains1,d2 in domains2]
end

splitdomain_sing(sing::SingularDiagonal, domain) =
    splitdomain_diagonal(domain)

splitdomain_diagonal(domain::ProductDomain) = splitdomain_diagonal(domain, components(domain)...)

function splitdomain_diagonal(domain::ProductDomain, domain1::AbstractInterval, domain2::AbstractInterval)
    T = prectype(domain)
    diff1 = domain1 \ domain2
    diff2 = domain2 \ domain1
    overlap = domain1 ∩ domain2
    if width(overlap) < 100eps(T)
        # "if isempty(overlap)" is not sufficient, because the overlap may be numerically small
        [domain1 × domain2]
    else
        a, b = extrema(overlap)
        if isempty(diff1) && isempty(diff2)
            [LowerRightTriangle(a, b), UpperRightTriangle(a, b)]
        elseif isempty(diff1)
            [LowerRightTriangle(a, b), UpperRightTriangle(a, b), domain1 × diff2]
        elseif isempty(diff2)
            [diff1 × domain2, LowerRightTriangle(a, b), UpperRightTriangle(a, b)]
        else
            [diff1 × domain2, LowerRightTriangle(a, b), UpperRightTriangle(a, b), overlap × diff2]
        end
    end
end


text:
using PSCOPF

using Test
using Dates
using DataStructures

@testset verbose=true  "test_seq_generator" begin

    #=
    INPUT :
        Network := Electric grid description
        TS := target time steps (dates d'intérêt)
        Mode de gestion = mode 1
        ECH := list of horizon points as generated by the generate_ech responsibility
    EXPECTED OUTPUT :
        At each horizon point, a sequence (an ordered list of operations i.e. AbstractRunnable)
    CARE POINTs:
        - generate_sequences assumes TS is sorted and unique ?
        - generate_sequences assumes ECH is sorted and unique ?
    =#

    @testset "mode_1" begin
        network = PSCOPF.Networks.Network()
        ts1 = Dates.DateTime("2015-01-01T11:00:00")
        TS = PSCOPF.create_target_timepoints(ts1)
        mode = PSCOPF.PSCOPF_MODE_1
        ECH = PSCOPF.generate_ech(network, TS, mode)

        sequence = PSCOPF.generate_sequence(network, TS, ECH, mode)

        EXPECTED_OPERATIONS = SortedDict(
            ts1 - Dates.Hour(4)      => [PSCOPF.EnergyMarket, PSCOPF.TSOOutFO],
            ts1 - Dates.Hour(1)      => [PSCOPF.EnergyMarketAtFO, PSCOPF.EnterFO, PSCOPF.TSOInFO],
            ts1 - Dates.Minute(30)   => [PSCOPF.TSOInFO],
            ts1 - Dates.Minute(15)   => [PSCOPF.TSOInFO],
            ts1                      => [PSCOPF.Assessment]
        )

        @test length(sequence) == length(ECH)
        for (ech, steps) in sequence.operations
            for (index,step) in enumerate(steps)
                @test isa(step, EXPECTED_OPERATIONS[ech][index])
            end
        end
    end

    @testset "mode_2" begin
        network = PSCOPF.Networks.Network()
        ts1 = Dates.DateTime("2015-01-01T11:00:00")
        TS = PSCOPF.create_target_timepoints(ts1)
        mode = PSCOPF.PSCOPF_MODE_2
        ECH = PSCOPF.generate_ech(network, TS, mode)

        sequence = PSCOPF.generate_sequence(network, TS, ECH, mode)

        EXPECTED_OPERATIONS = SortedDict(
            ts1 - Dates.Hour(4)      => [PSCOPF.EnergyMarket, PSCOPF.TSOOutFO],
            ts1 - Dates.Hour(1)      => [PSCOPF.EnergyMarketAtFO, PSCOPF.EnterFO, PSCOPF.TSOOutFO],
            ts1 - Dates.Minute(30)   => [PSCOPF.EnergyMarket, PSCOPF.TSOBilevel, PSCOPF.BalanceMarket],
            ts1 - Dates.Minute(15)   => [PSCOPF.EnergyMarket, PSCOPF.TSOBilevel, PSCOPF.BalanceMarket],
            ts1                      => [PSCOPF.Assessment]
        )

        @test length(sequence) == length(ECH)
        for (ech, steps) in sequence.operations
            for (index,step) in enumerate(steps)
                @test isa(step, EXPECTED_OPERATIONS[ech][index])
            end
        end
    end

    @testset "mode_3" begin
        network = PSCOPF.Networks.Network()
        ts1 = Dates.DateTime("2015-01-01T11:00:00")
        TS = PSCOPF.create_target_timepoints(ts1)
        mode = PSCOPF.PSCOPF_MODE_3
        ECH = PSCOPF.generate_ech(network, TS, mode)

        sequence = PSCOPF.generate_sequence(network, TS, ECH, mode)

        EXPECTED_OPERATIONS = SortedDict(
            ts1 - Dates.Hour(4)      => [PSCOPF.EnergyMarket, PSCOPF.TSOOutFO],
            ts1 - Dates.Hour(1)      => [PSCOPF.EnergyMarket, PSCOPF.EnterFO, PSCOPF.TSOBilevel], #FIXME or maybe TSOBilevelAtFO
            ts1 - Dates.Minute(30)   => [PSCOPF.BalanceMarket],
            ts1 - Dates.Minute(15)   => [PSCOPF.BalanceMarket],
            ts1                      => [PSCOPF.Assessment]
        )

        @test length(sequence) == length(ECH)
        for (ech, steps) in sequence.operations
            for (index,step) in enumerate(steps)
                @test isa(step, EXPECTED_OPERATIONS[ech][index])
            end
        end
    end

end

max_stars_repo_name: klorel/pscopf

max_issues_repo_name: klorel/pscopf

max_forks_repo_name: klorel/pscopf


text:
using Flux, Base.Test

srand(0)

@testset "Flux" begin

include("utils.jl")
include("tracker.jl")
include("layers/normalisation.jl")
include("layers/stateless.jl")
include("optimise.jl")
include("data.jl")

if Base.find_in_path("CuArrays") ≠ nothing
  include("cuda/cuda.jl")
end

end


text:
module CompressedSensing

using LinearAlgebra
using SparseArrays

using LazyInverses
using WoodburyFactorizations
using UpdatableQRFactorizations

const AbstractMatOrFac{T} = Union{AbstractMatrix{T}, Factorization{T}}

using Base.Threads
using StatsBase: sample, mean

using JuMP
using JuMP: Model, optimizer_with_attributes, @variable, @objective, @constraint
using JuMP: SecondOrderCone
using Clp # for linear program
using ECOS # ECOS is interior point method (more accurate but slower on large problems)
# using SCS # operator splitting method, fast but not high accuracy

abstract type Update{T} end
(U::Update)(x) = update!(U, x)

# Compressed Sensing algorithms
include("util.jl") # sparse data generator
include("oblivious.jl")
include("matchingpursuit.jl")
include("forward.jl")
include("backward.jl")
include("twostage.jl")
include("stepwise.jl")
include("sbl.jl")
include("basispursuit.jl")

end


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

include("model_CloudError.jl")
include("model_CloudErrorBody.jl")
include("model_DnsNameAvailabilityResult.jl")


text:
import Mustache

#Default options
const defaultParams =
      Dict{Symbol,Any}(:storeresults => false,
                                :doc_number => 0,
                                :chunk_defaults =>
                                Dict{Symbol,Any}(
                                :echo=> true,
                                :results=> "markup",
                                :hold => false,
                                :fig=> true,
                                :include=> true,
                                :eval => true,
                                :tangle => true,
                                :cache => false,
                                :fig_cap=> nothing,
                                #Size in inches
                                :fig_width => 6,
                                :fig_height => 4,
                                :fig_path=> "figures",
                                :dpi => 96,
                                :term=> false,
                                :display => false,
                                :prompt => "\njulia> ",
                                :label=> nothing,
                                :wrap=> true,
                                :line_width => 75,
                                :engine=> "julia",
                                #:option_AbstractString=> "",
                                #Defined in formats
                                :fig_ext => nothing,
                                :fig_pos=> nothing,
                                :fig_env=> nothing,
                                :out_width=> nothing,
                                :out_height=> nothing,
                                :skip=>false
                                )
                            )
#This one can be changed at runtime, initially a copy of defaults
const rcParams = deepcopy(defaultParams)

"""
    set_chunk_defaults(opts::Dict{Symbol, Any})

Set default options for code chunks, use [`get_chunk_defaults`](@ref) to see the current values.

E.g.: set default `dpi` to 200 and `fig_width` to 8

```julia
julia> set_chunk_defaults(Dict{Symbol, Any}(:dpi => 200, fig_width => 8))
```
"""
function set_chunk_defaults(opts::Dict{Symbol, Any})
  merge!(rcParams[:chunk_defaults], opts)
  return nothing
end

"""
`get_chunk_defaults()`

Get default options used for code chunks.
"""
function get_chunk_defaults()
  return(rcParams[:chunk_defaults])
end

"""
`restore_chunk_defaults()`

Restore Weave.jl default chunk options
"""
function restore_chunk_defaults()
  rcParams[:chunk_defaults] = defaultParams[:chunk_defaults]
  return nothing
end

"""Combine format specific and common options from document header"""
function combine_args(args, doctype)
    common = Dict()
    specific = Dict()
    for key in keys(args)
        if key ∈ keys(Weave.formats)
            specific[key] = args[key]
        else
            common[key] = args[key]
        end
    end
    haskey(specific, doctype) && merge!(common, specific[doctype])
    common
end

getvalue(d::Dict, key , default) = haskey(d, key) ? d[key] : default

"""
header_args(doc::WeaveDoc)`

Get weave arguments from document header
"""
function header_args(doc::WeaveDoc, out_path, mod, fig_ext, fig_path,
                            cache_path, cache, throw_errors,template,
                            highlight_theme, css,
                            pandoc_options, latex_cmd)
    args = getvalue(doc.header, "options", Dict())
    doctype = getvalue(args, "doctype", doc.doctype)
    args = combine_args(args, doctype)
    informat = getvalue(args, "informat", :auto)
    out_path = getvalue(args, "out_path", out_path)
    out_path == ":pwd" && (out_path = :pwd)
    isa(out_path, Symbol) || (out_path = joinpath(dirname(doc.source), out_path))
    mod = Symbol(getvalue(args, "mod", mod))
    fig_path = getvalue(args, "fig_path", fig_path)
    fig_ext = getvalue(args, "fig_ext", fig_ext)
    cache_path = getvalue(args, "cache_path", cache_path)
    cache = Symbol(getvalue(args, "cache", cache))
    throw_errors = getvalue(args, "throw_errors", throw_errors)
    template = getvalue(args, "template", template)
    if template != nothing && !isa(template, Mustache.MustacheTokens) && !isempty(template)
        template = joinpath(dirname(doc.source), template)
    end
    highlight_theme = getvalue(args, "highlight_theme", highlight_theme)
    css = getvalue(args, "css", css)
    pandoc_options = getvalue(args, "pandoc_options", pandoc_options)
    latex_cmd = getvalue(args, "latex_cmd", latex_cmd)

    return (doctype, informat, out_path, args, mod, fig_path, fig_ext,
      cache_path, cache, throw_errors, template, highlight_theme, css,
      pandoc_options, latex_cmd)
end

"""
`header_chunk_defaults!(doc::WeaveDoc)`

Get chunk defaults from header and update
"""
function header_chunk_defaults!(doc::WeaveDoc)
    for key in keys(doc.chunk_defaults)
        if haskey(doc.header["options"], String(key))
             doc.chunk_defaults[key] = doc.header["options"][String(key)]
        end
    end
end


text:
"""
    xnew,unew,cnew,sigmanew = forward_pass(traj_new, x0,u,x,α,f,costfun,lims,diff)
# Arguments
- α: step size (αk is applied to old trajectory)
- diff: function to determine difference `diff(xnew[:,i], x[:,i])`
- f: forward dynamics `x(k+1)  = f(x(k), u(k), k)`
- `cnew = costfun(xnew, unew)`
"""
@views function forward_pass(traj_new, x0,u,x,α,f,costfun,lims,diff)
    n         = size(x0,1)
    m,N       = size(u)
    xnew      = Array{eltype(x0)}(undef,n,N)
    xnew[:,1] = x0
    unew      = copy(u)
    cnew      = zeros(N)
    for i = 1:N
        if !isempty(traj_new)
            unew[:,i] .+= traj_new.k[:,i]*α
            dx = diff(xnew[:,i], x[:,i])
            unew[:,i] .+= traj_new.K[:,:,i]*dx
        end
        if !isempty(lims)
            unew[:,i] = clamp.(unew[:,i],lims[:,1], lims[:,2])
        end
        xnewi  = f(xnew[:,i], unew[:,i], i)
        if i < N
            xnew[:,i+1] = xnewi
        end
    end
    cnew = costfun(xnew, unew)

    return xnew,unew,cnew
end

@views function forward_covariance(model, x, u, traj)
    fx,fu,fxx,fxu,fuu = df(model, x, u)
    n        = size(fx,1)
    m        = size(fu,2)
    N        = size(fx,3)
    R1       = covariance(model,x,u) # Simple empirical prediction covariance
    Σ0       = R1 # TODO: I was lazy here
    ix = 1:n
    iu = n+1:n+m
    sigmanew = Array{eltype(x)}(undef,n+m,n+m,N)
    sigmanew[ix,ix,1] = Σ0
    for i = 1:N-1
        K,Σ = traj.K[:,:,i], traj.Σ[:,:,i]
        sigmanew[ix,ix,i+1] = fx[:,:,i]*sigmanew[ix,ix,i]*fx[:,:,i]' + R1 # Iterate dLyap forward
        sigmanew[iu,ix,i] = K*sigmanew[ix,ix,i]
        sigmanew[ix,iu,i] = sigmanew[ix,ix,i]*K'
        sigmanew[iu,iu,i] = K*sigmanew[ix,ix,i]*K' + Σ
    end
    sigmanew
end


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

@doc raw"""FSGroupStrategyOptions defines the strategy type and options used to create the strategy.

    IoK8sApiPolicyV1beta1FSGroupStrategyOptions(;
        ranges=nothing,
        rule=nothing,
    )

    - ranges::Vector{IoK8sApiPolicyV1beta1IDRange} : ranges are the allowed ranges of fs groups.  If you would like to force a single fs group then supply a single range with the same start and end. Required for MustRunAs.
    - rule::String : rule is the strategy that will dictate what FSGroup is used in the SecurityContext.
"""
mutable struct IoK8sApiPolicyV1beta1FSGroupStrategyOptions <: SwaggerModel
    ranges::Any # spec type: Union{ Nothing, Vector{IoK8sApiPolicyV1beta1IDRange} } # spec name: ranges
    rule::Any # spec type: Union{ Nothing, String } # spec name: rule

    function IoK8sApiPolicyV1beta1FSGroupStrategyOptions(;ranges=nothing, rule=nothing)
        o = new()
        validate_property(IoK8sApiPolicyV1beta1FSGroupStrategyOptions, Symbol("ranges"), ranges)
        setfield!(o, Symbol("ranges"), ranges)
        validate_property(IoK8sApiPolicyV1beta1FSGroupStrategyOptions, Symbol("rule"), rule)
        setfield!(o, Symbol("rule"), rule)
        o
    end
end # type IoK8sApiPolicyV1beta1FSGroupStrategyOptions

const _property_map_IoK8sApiPolicyV1beta1FSGroupStrategyOptions = Dict{Symbol,Symbol}(Symbol("ranges")=>Symbol("ranges"), Symbol("rule")=>Symbol("rule"))
const _property_types_IoK8sApiPolicyV1beta1FSGroupStrategyOptions = Dict{Symbol,String}(Symbol("ranges")=>"Vector{IoK8sApiPolicyV1beta1IDRange}", Symbol("rule")=>"String")
Base.propertynames(::Type{ IoK8sApiPolicyV1beta1FSGroupStrategyOptions }) = collect(keys(_property_map_IoK8sApiPolicyV1beta1FSGroupStrategyOptions))
Swagger.property_type(::Type{ IoK8sApiPolicyV1beta1FSGroupStrategyOptions }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiPolicyV1beta1FSGroupStrategyOptions[name]))}
Swagger.field_name(::Type{ IoK8sApiPolicyV1beta1FSGroupStrategyOptions }, property_name::Symbol) =  _property_map_IoK8sApiPolicyV1beta1FSGroupStrategyOptions[property_name]

function check_required(o::IoK8sApiPolicyV1beta1FSGroupStrategyOptions)
    true
end

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


text:
# Approximations
include("./approximations/gradient_descent.jl")
include("./approximations/SR1.jl")
include("./approximations/DFP.jl")
include("./approximations/BB.jl")
include("./approximations/BFGS.jl")
include("./approximations/DBFGS.jl")
include("./approximations/LBFGS.jl")
include("./approximations/newton.jl")
include("./approximations/scalings.jl")

include("update_qn.jl")
include("find_direction.jl")


text:
"""
# script test_illusions

- Julia version: 1.4
- Author: niallcullinane
- Date: 2020-08-20

Script to run set of images with LAMINART and plot layers and activation.
"""

using DrWatson
@quickactivate "MEngProject"
using MEngProject,
    CUDA,
    DifferentialEquations,
    PyPlot,
    NNlib,
    ImageFiltering,
    Images,
    ImageIO,
    MEngProject.LaminartKernels,
    MEngProject.LaminartInitFunc,
    MEngProject.Utils,
    BenchmarkTools,
    Test

using OrdinaryDiffEq,
    ParameterizedFunctions, LSODA, Sundials, DiffEqDevTools, Noise

batch = 9003

files = [
    "Iine_gap_1_100_gs.png",
    "Iine_gap_2_100_gs.png",
    "Iine_gap_3_100_gs.png",
    "Iine_gap_4_100_gs.png",
    "Iines_gaps_100_gs.png",
    "diag_dots_100_gs.png",
    "diag_gap_100_gs.png",
    "kan_sq_cont.png",
    "kan_sq_cont_l.png",
    "mo05709.png",
    "stairs_100gs.png",
    "stairs_200gs.png",
    "viper00187.png",
    "viper00661.jpg",
    "viper00715.png",
    "viper00717.png",
    "viper00720.png",
    "viper00721.png",
    "viper00842.jpg",
    "viper00891.jpg",
    "viper00904.jpg",
    "viper00921.jpg",
    "viper01006.jpg",
    "viper01333.jpg",
]

batch_ = string(batch, "_", rand(1000:9999));
mkdir(plotsdir(string("input_plot", batch_)));

tspan = (0.0f0, 800f0);

for file in files
    p = LaminartInitFunc.parameterInit_conv_cpu(
        datadir("img", file),
        Parameters.parameters_f32,
    )

    Utils.plot_gs(
        p.I[:, :, 1, 1],
        name = file,
        loc = plotsdir(string("input_plot", batch_)),
        clbar = true,
        axMax = 1,
        save = true,
    )
end


text:
# nvprof --profile-from-start off julia /path/to/ncu.jl
# or
# ncu --mode=launch julia /path/to/ncu.jl
using GPUAcceleratedTracking, CUDA, Tracking, GNSSSignals, StructArrays;
import Tracking: Hz, ms;
enable_gpu = Val(true)
num_samples = 50000
sampling_frequency = num_samples / 1ms |> Hz
num_ants = 1
num_correlators = 3
system = GPSL1(use_gpu = enable_gpu)
codes = system.codes
#convert to text_mem
codes = CuTexture(
    CuTextureArray(codes),
    address_mode = CUDA.ADDRESS_MODE_WRAP,
    interpolation = CUDA.NearestNeighbour(),
    normalized_coordinates = true
)
code_frequency = get_code_frequency(system)
code_length = get_code_length(system)
start_code_phase = 0.0f0
prn = 1
correlator = EarlyPromptLateCorrelator(NumAnts(num_ants), NumAccumulators(num_correlators))
correlator_sample_shifts = get_correlator_sample_shifts(system, correlator, sampling_frequency, 0.5)
num_of_shifts = correlator_sample_shifts[end] - correlator_sample_shifts[1]
code_replica = CUDA.zeros(Float32, num_samples + num_of_shifts)
code_replica_kernel = @cuda launch=false gen_code_replica_texture_mem_strided_kernel!(
            code_replica,
            codes,
            code_frequency,
            sampling_frequency,
            start_code_phase,
            prn,
            num_samples,
            num_of_shifts,
            code_length
        )
blocks_per_grid, threads_per_block = launch_configuration(code_replica_kernel.fun)
CUDA.@profile @cuda threads=threads_per_block blocks=blocks_per_grid gen_code_replica_texture_mem_strided_kernel!(
    code_replica,
    codes,
    code_frequency,
    sampling_frequency,
    start_code_phase,
    prn,
    num_samples,
    num_of_shifts,
    code_length
)

max_stars_repo_name: mfkiwl/GPUAcceleratedTracking

max_issues_repo_name: mfkiwl/GPUAcceleratedTracking

max_forks_repo_name: mfkiwl/GPUAcceleratedTracking


text:
import AbstractPlotting.convert_arguments

# bin centers
function convert_arguments(P::Type{<:Scatter}, B::SSHBinner)
    convert_arguments(P, bin_positions(B))
end
# mesh with vertices @ bin centers (makes color mapping easy)
function convert_arguments(P::Type{<:AbstractPlotting.Mesh}, B::SSHBinner)
    convert_arguments(P, vertex_mesh(B))
end
# boundaries of bins (slightly enlarged to work with mesh)
function convert_arguments(P::Type{<:LineSegments}, B::SSHBinner)
    convert_arguments(P, 1.001line_segments(B))
end

# Not sure if I wanna keep this. Would probably be nice to have a colorlegend
# by default here
"""
    histogram(B::SSHBinner)
    histogram!(scene, B::SSHBinner)

Plots a histogram for the given SSHBinner.

# Attributes

- `outline = false`: Draws lines around each bin
- `colormap = :matter`: Colormap used to visualize bins
- `colorscale = :absolute`: Scale used when mapping bins to colors. Can be
:relative (sets colorrange = extrema(B.bins)), :absolute (sets only the max) or
:static (colorrange remains unchanged)
- `colorrange = [0, maxiumum(b.bins)]`: Range of bin fillings represented by the
colormap.
- `transparency = false`: Set this to true if your colormap is tranparent
- `linecolor = :black`: Color of outline.
"""
@recipe(Histogram, B) do scene
    Theme(
        outline = false,
        colormap = :matter,
        transparency = false,
        linecolor = :black,
        colorscale = :absolute,
        colorrange = (0, 1)
    )
end

function map_bin_to_color(binner, colormap, colorscale, colorrange)
    map(binner, colormap, colorscale) do binner, cm, scale
        ex = extrema(binner.bins)
        if scale == :absolute
            colorrange[] = (0, ex[2])
        elseif scale == :static
        else
            colorrange[] = (ex[1], ex[2])
        end
        map(binner.bins) do b
            AbstractPlotting.RGBAf0(AbstractPlotting.AbstractPlotting.interpolated_getindex(
                to_colormap(cm), Float64(b), colorrange[]
            ))
        end
    end
end

function AbstractPlotting.plot!(plot::Histogram)
    colors = map_bin_to_color(
        plot[:B], plot[:colormap], plot[:colorscale], plot[:colorrange]
    )

    mesh!(plot, plot[:B], color = colors, transparency = plot[:transparency])
    linesegments!(plot, plot[:B], color = plot[:linecolor], visible=plot[:outline])

    plot
end

export histogram, histogram!

# """
#     histogram(B::SSHBinner[, color_method = to_hue])
#
# Plots a sphere where each bin filling is mapped to a color through color_method.
# """
# function histogram(B::SSHBinner, color_method=to_hue)
#     _colors = color_method(B)
#     is_transparent = any(c -> c.alpha < 1.0, _colors)
#     AbstractPlotting.mesh(vertex_mesh(B), color = _colors, transparency = is_transparent)
# end
#
#
# """
#     plot_debug(B::SSHBinner[, color_method = to_hue])
#
# Plots a sphere with colors representing bin fillings, approximate bin area
# (rectangles) and points at the bin centers.
# """
# function plot_debug(B::SSHBinner, color_method = to_hue)
    # _colors = color_method(B)
    # is_transparent = any(c -> c.alpha < 1.0, _colors)
    # scene = AbstractPlotting.mesh(vertex_mesh(B), color = _colors, transparency = is_transparent)
    # AbstractPlotting.linesegments!(scene, 1.001line_segments(B), color = :black)
    # # AbstractPlotting.mesh!(scene1, to_rects(B, extrude), color = RGBA(1, 1, 1, 0.2), transparency = true)
    # AbstractPlotting.scatter!(scene, 1.005bin_positions(B), color=:black, markersize=0.01)
    # scene
# end


text:
# Autogenerated wrapper script for libvorbis_jll for x86_64-apple-darwin
export libvorbis, libvorbisenc, libvorbisfile

using Ogg_jll
JLLWrappers.@generate_wrapper_header("libvorbis")
JLLWrappers.@declare_library_product(libvorbis, "@rpath/libvorbis.0.dylib")
JLLWrappers.@declare_library_product(libvorbisenc, "@rpath/libvorbisenc.2.dylib")
JLLWrappers.@declare_library_product(libvorbisfile, "@rpath/libvorbisfile.3.dylib")
function __init__()
    JLLWrappers.@generate_init_header(Ogg_jll)
    JLLWrappers.@init_library_product(
        libvorbis,
        "lib/libvorbis.0.dylib",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@init_library_product(
        libvorbisenc,
        "lib/libvorbisenc.2.dylib",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@init_library_product(
        libvorbisfile,
        "lib/libvorbisfile.3.dylib",
        RTLD_LAZY | RTLD_DEEPBIND,
    )

    JLLWrappers.@generate_init_footer()
end  # __init__()


text:
@testset "wrapping_integers_wrap.jl" begin
    @test wrap(3 * (UInt(1) << 32), WrappingInt32(0)) == WrappingInt32(0)
    @test wrap(3 * (UInt(1) << 32) + 17, WrappingInt32(15)) == WrappingInt32(32)
    @test wrap(7 * (UInt(1) << 32) - 2, WrappingInt32(15)) == WrappingInt32(13)
end

max_stars_repo_name: AquaIndigo/JLSponge

max_issues_repo_name: AquaIndigo/JLSponge

max_forks_repo_name: AquaIndigo/JLSponge


text:
export ChunkedVector

immutable ChunkedVector{T,A<:AVector,N} <: AVector{T}
  data::Vector{A}
end

typealias ChunkedVectorT{T,A,N} Type{ChunkedVector{T,A,N}}

(::ChunkedVectorT{T,A,N}){T,A,N}(data = A[A()]) =
  ChunkedVector{T,A,N}(data)

(::ChunkedVectorT{T,A}){T,A}(data = A[A()]) =
  ChunkedVector{T,A,128}(data)

(::ChunkedVectorT{T}){T}(a...) =
  ChunkedVector{T,Vector{T}}(a...)

chunksize{T,A,N}(::Type{ChunkedVector{T,A,N}}) = N
chunksize(xs::ChunkedVector) = chunksize(typeof(xs))

chunks(xs::ChunkedVector) = xs.data

function makeroom!(xs::ChunkedVector)
  push!(xs.data, similar(xs.data[end], 0))
end

function Base.push!(xs::ChunkedVector, y)
  length(xs.data[end]) == chunksize(xs) && makeroom!(xs)
  push!(xs.data[end], y)
  return xs
end

function ChunkedVector(xs, a...)
  v = ChunkedVector{eltype(xs)}(a...)
  for x in xs
    push!(v, x)
  end
  return v
end

function index_slow(xs::ChunkedVector, i)
  j, i′ = divrem(i-1, chunksize(xs))
  j + 1, i′ + 1
end

@generated function index_fast(xs::ChunkedVector, i)
  pow = round(Int, log2(chunksize(xs)))
  :((i-1) >> $pow + 1, (i-1) & $(2^pow-1) + 1)
end

@generated function index(xs::ChunkedVector, i)
  ispow2(chunksize(xs)) ?
    :(index_fast(xs, i)) :
    :(index_slow(xs, i))
end

@inline function Base.getindex(xs::ChunkedVector, i::Integer)
  j, i′ = index(xs, i)
  @inbounds return xs.data[j][i′]
end

function Base.setindex!(xs::ChunkedVector, v, i::Integer)
  j, i′ = index(xs, i)
  @inbounds xs.data[j][i′] = v
end

Base.length(xs::ChunkedVector) = chunksize(xs)*(length(xs.data)-1)+length(xs.data[end])

Base.size(xs::ChunkedVector) = (length(xs),)

# Iteration utils

immutable ChunkIter{N}
  max::Tuple{Int,Int}
end

chunksize{N}(::Type{ChunkIter{N}}) = N
chunksize(it::ChunkIter) = chunksize(typeof(it))

Base.start(::ChunkIter) = (1, 0)

Base.done(it::ChunkIter, state) =
  state[2] == it.max[2] && state[1] == it.max[1]

@inline function _next(it::ChunkIter, last)
  i, j = last
  j == chunksize(it) ? (i+1, 1) : (i, j+1)
end

function Base.next(it::ChunkIter, last)
  next = _next(it, last)
  next, next
end

Base.length(it::ChunkIter) = (it.max[1]-1)*chunksize(it)+it.max[2]

# ChunkVector iteration

ChunkIter(xs::ChunkedVector) =
  ChunkIter{chunksize(xs)}((length(xs.data), length(xs.data[end])))

@iter xs::ChunkedVector -> ChunkIter(xs)

@inline function Base.next(xs::ChunkedVector, sub::SubIter)
  i, j = _next(sub.iter, sub.state)
  @inbounds return xs.data[i][j], SubIter(sub.iter, (i, j))
end