text:
# ---
# title: 36. Valid Sudoku
# id: problem36
# author: Tian Jun
# date: 2020-10-31
# difficulty: Medium
# categories: Hash Table
# link: <https://leetcode.com/problems/valid-sudoku/description/>
# hidden: true
# ---
#
# Determine if a `9 x 9` Sudoku board is valid. Only the filled cells need to be
# validated  **according to the following rules** :
#
#   1. Each row must contain the digits `1-9` without repetition.
#   2. Each column must contain the digits `1-9` without repetition.
#   3. Each of the nine `3 x 3` sub-boxes of the grid must contain the digits `1-9` without repetition.
#
# **Note:**
#
#   * A Sudoku board (partially filled) could be valid but is not necessarily solvable.
#   * Only the filled cells need to be validated according to the mentioned rules.
#
#
#
# **Example 1:**
#
# ![](https://upload.wikimedia.org/wikipedia/commons/thumb/f/ff/Sudoku-
# by-L2G-20050714.svg/250px-Sudoku-by-L2G-20050714.svg.png)
#
#
#
#     Input: board =
#     [["5","3",".",".","7",".",".",".","."]
#     ,["6",".",".","1","9","5",".",".","."]
#     ,[".","9","8",".",".",".",".","6","."]
#     ,["8",".",".",".","6",".",".",".","3"]
#     ,["4",".",".","8",".","3",".",".","1"]
#     ,["7",".",".",".","2",".",".",".","6"]
#     ,[".","6",".",".",".",".","2","8","."]
#     ,[".",".",".","4","1","9",".",".","5"]
#     ,[".",".",".",".","8",".",".","7","9"]]
#     Output: true
#
#
# **Example 2:**
#
#
#
#     Input: board =
#     [["8","3",".",".","7",".",".",".","."]
#     ,["6",".",".","1","9","5",".",".","."]
#     ,[".","9","8",".",".",".",".","6","."]
#     ,["8",".",".",".","6",".",".",".","3"]
#     ,["4",".",".","8",".","3",".",".","1"]
#     ,["7",".",".",".","2",".",".",".","6"]
#     ,[".","6",".",".",".",".","2","8","."]
#     ,[".",".",".","4","1","9",".",".","5"]
#     ,[".",".",".",".","8",".",".","7","9"]]
#     Output: false
#     Explanation: Same as Example 1, except with the **5** in the top left corner being modified to **8**. Since there are two 8's in the top left 3x3 sub-box, it is invalid.
#
#
#
#
# **Constraints:**
#
#   * `board.length == 9`
#   * `board[i].length == 9`
#   * `board[i][j]` is a digit or `'.'`.
#
#
## @lc code=start
using LeetCode

## add your code here:
## @lc code=end


text:
# ------------ includes ------------ #

include("foreign_model_interfaces/blackbox.jl")
include("foreign_model_interfaces/soss.jl")
include("foreign_model_interfaces/turing.jl")
include("foreign_model_interfaces/gen.jl")
include("foreign_model_interfaces/flux.jl")

function primitive end
function load_soss_fmi end
function load_turing_fmi end
function load_flux_fmi end
function load_gen_fmi end

foreign(addr::A, args...) where A <: Address = error("(foreign) call with address $addr evaluated outside of the tracer.\nThis normally occurs because you're not matching the dispatch correctly.")

# ------------ Documentation -------------- #

@doc(
"""
```julia
@primitive function logpdf(fn::typeof(foo), args, foo_ret)
    ...
end
```

`@primitive` is a convenience metaprogramming construct which derives contextual dispatch definitions for functions which you'd like the tracer to "summarize" to a single site. This is one mechanism which gives the user more control over the structure of the choice map structure in their programs.

Example:

```julia
geo(p::Float64) = rand(:flip, Bernoulli(p)) ? 1 : 1 + rand(:geo, geo, p)

# Define as primitive.
@primitive function logpdf(fn::typeof(geo), p, count)
    return Distributions.logpdf(Geometric(p), count)
end

function foo()
    ret = rand(:geo, geo, 0.3)
    ret
end
```

The above program would summarize the trace into a single choice site for `:geo`, as if `geo` was a primitive distribution.

```julia
  __________________________________

               Addresses

 geo : 42
  __________________________________
```
""", primitive)

@doc(
"""
```julia
foreign(addr::A, m, args...) where A <: Address
```

Activate a foreign model interface. The tracer will treat this is a specialized call site, depending on the type of `m`. Currently supports `typeof(m) <: Soss.Model` and `typeof(m) <: Gen.GenerativeFunction`.

""", foreign)

@doc(
"""
```julia
Jaynes.@load_soss_fmi
```

`@load_soss_fmi` loads the [Soss.jl](https://github.com/cscherrer/Soss.jl) foreign model interface extension. This allows you to utilize `Soss` models in your modeling.

Example:

```julia
Jaynes.@load_soss_fmi()

# A Soss model.
m = @model σ begin
    μ ~ Normal()
    y ~ Normal(μ, σ) |> iid(5)
end

bar = () -> begin
    x = rand(:x, Normal(5.0, 1.0))
    soss_ret = foreign(:foo, m, (σ = x,))
    return soss_ret
end
```

This interface currently supports all the inference interfaces (e.g. `simulate`, `generate`, `score`, `regenerate`, `update`, `propose`) which means that you can use any of the inference algorithms in the standard inference library.
""", load_soss_fmi)

@doc(
"""
```julia
Jaynes.@load_gen_fmi
```

`@load_gen_fmi` loads the [Gen.jl](https://www.gen.dev/) foreign model interface extension. This allows you to utilize `Gen` models (in any of Gen's DSLs) in your modeling.

Example:

```julia
Jaynes.@load_gen_fmi()

@gen (static) function foo(z::Float64)
    x = @trace(normal(z, 1.0), :x)
    y = @trace(normal(x, 1.0), :y)
    return x
end

Gen.load_generated_functions()

bar = () -> begin
    x = rand(:x, Normal(0.0, 1.0))
    return foreign(:foo, foo, x)
end

ret, cl = Jaynes.simulate(bar)
```

This interface currently supports all the inference interfaces (e.g. `simulate`, `generate`, `score`, `regenerate`, `update`, `propose`) which means that you can use any of the inference algorithms in the standard inference library.
""", load_gen_fmi)

ext: jl


text:
struct BreastCancer <: Name end

# mandatory methods
format(::Type{BreastCancer}) = TabularData
task(::Type{BreastCancer}) = TwoClass
source(::Type{BreastCancer}) = "https://archive.ics.uci.edu/ml/datasets/Breast+Cancer"

function downloadlink(::Type{BreastCancer})
    return "https://archive.ics.uci.edu/ml/machine-learning-databases/breast-cancer/breast-cancer.data"
end

function preprocess(N::Type{BreastCancer})
    return path -> preprocess_csv(N, path, :train; istarget = 1, iscategorical = 2:10)
end

nattributes(::Type{BreastCancer}) = (9, )
ninstances(::Type{BreastCancer}) = (286, 0, 0)
positivelabel(::Type{BreastCancer}) = "recurrence-events"

# optional methods
function checksum(::Type{BreastCancer})
    return "ca7d3fa97b62ff967b6894ffbb3acaaf4a51062e0149e5950b3ad6f685970b65"
end

name(::Type{BreastCancer}) = "Breast Cancer Data Set"
function author(::Type{BreastCancer})
    return ["Matjaz Zwitter", "Milan Soklic", "Ming Tan", "Jeff Schlimmer"]
end

max_stars_repo_name: VaclavMacha/DatasetProvider

max_issues_repo_name: VaclavMacha/DatasetProvider

max_forks_repo_name: VaclavMacha/DatasetProvider


text:
using CutCellDG
using Test
include("../useful_routines.jl")

x0 = [-1.0]
meshwidth = [2.0]
nelements = [3]
nodesperelement = 4

cgmesh = CutCellDG.CGMesh(x0,meshwidth,nelements,nodesperelement)
nc = CutCellDG.nodal_connectivity(cgmesh)
testnc = reshape([1,2,3,4,4,5,6,7,7,8,9,10],4,:)

@test allequal(nc,testnc)


text:
using Parameters, Plots, VAWTools
include("../itmix-setup.jl")

if !isdefined(:gl_all_radar)
    glacier = :Unteraar
    gid = BM.ITMIXGlacier(glacier, 2)
    gl_all_radar,gb,pp,pm,pn,pl = BM.init_forward(gid)
end

pm = BM.MPara(pm,
              sigma_h_model=20,
              error_model=[BM.ErrorModel.sqrerr, BM.ErrorModel.abserr][1],
              nthin_h=1)

runtyp = [:test,:testmid,:prodlow,:prodmid,:prod][3]
fit_target = [BM.FitTarget.h, BM.FitTarget.iv, BM.FitTarget.h_iv, BM.FitTarget.length][1]
n_1d_btilde = 3
n_1d_fsl = 3
n_1d_temp = 1
rm_iv = false
r_tracks = [5,34,42]

(gln, gbn, th0d,
 theta0, logposterior, logposterior1d,
 logprior, logprior_debug, loglikelihood, loglikelihood1d,
 fwdm_fn, fwdm1d_fn, pmcmc_defaults, fit_target) =
     setup_one_exp(r_tracks, gl_all_radar, gb, pp, pm, pn,
                   n_1d_btilde, n_1d_fsl, n_1d_temp,
                   runtyp, fit_target, rm_iv)

theta = zeros(theta0.th0)
theta[1:3] = [-2.5, 0.7, 0.0]

xtheta = VAWTools.piecewiselinear(gb.ele, gb.xmid).(theta0.th_ele[:btilde])

new_btilde = BM.make_fields1d(theta, theta0)[1]

plot(gb.xmid/1e3, gb.bdot-gb.dhdt, xlabel="χ (km)", ylabel="b̃ (ma⁻¹)", label="Measurements", size=(400,300))
plot!(xtheta/1e3, theta[1:3], marker=:circle, label="θ")
plot!(gb.xmid/1e3, new_btilde, label="Used")

#savefig("figs/explanation-1D-vars.png")
#savefig("figs/explanation-1D-vars.pdf")


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_Usage.jl")
include("model_UsageName.jl")
include("model_UsagesListResult.jl")


text:
struct GAE{T,S}
    encoder::T
    decoder::S
end

GAE(enc, σ::Function=identity) = GAE(enc, InnerProductDecoder(σ))

@treelike GAE

function (g::GAE)(X::AbstractMatrix)
    Z = g.encoder(X)
    A = g.decoder(Z)
    A
end

struct VGAE{T,S}
    encoder::T
    decoder::S
end

function VGAE(enc, h_dim::Integer, z_dim::Integer, σ::Function=identity)
    VGAE(VariationalEncoder(enc, h_dim, z_dim), InnerProductDecoder(σ))
end

@treelike VGAE

function (g::VGAE)(X::AbstractMatrix)
    Z = g.encoder(X)
    A = g.decoder(Z)
    A
end

struct InnerProductDecoder
    σ
end

@treelike InnerProductDecoder

(i::InnerProductDecoder)(Z::AbstractArray) = i.σ(Z'*Z)

struct VariationalEncoder
    nn
    μ
    logσ
    z_dim::Integer

    function VariationalEncoder(nn, h_dim::Integer, z_dim::Integer)
        new(nn, Dense(h_dim, z_dim), Dense(h_dim, z_dim), z_dim)
    end
end

@treelike VariationalEncoder

function (ve::VariationalEncoder)(X::AbstractMatrix)
    h = ve.nn(X)
    μ, logσ = summarize(ve, h)
    Z = sampling(ve, μ, logσ)
    Z
end

summarize(ve::VariationalEncoder, h::AbstractMatrix) = (ve.μ(h), ve.logσ(h))
sampling(::VariationalEncoder, μ, logσ) = μ .+ exp.(logσ) .* randn(size(logσ))

# function logpdf(b::Bernoulli, y::Bool; T::Real=Float32)
#     y * log(b.p + eps(T)) + (one(T) - y) * log(one(T) - b.p + eps(T))
# end
#
# # KL-divergence between approximation posterior and N(0, 1) prior.
# kl_q_p(μ, logσ; T::Real=Float32) = T(0.5) * sum(exp.(T(2) .* logσ) + μ.^2 .- one(T) .+ logσ.^2)
#
# # logp(x|z) - decoder
# logp_x_z(x, z) = sum(logpdf.(Bernoulli.(f(z)), x))
#
# function loss(X; T=eltype(X), β=one(T), λ=T(0.01))
#     N = size(X, 1)  # batch size
#     μ̂, logσ̂ = g(X)
#     L̄ = 1//N * (logp_x_z(X, sampling.(μ̂, logσ̂)) - β * kl_q_p(μ̂, logσ̂))
#     -L̄ + λ * sum(x->sum(x.^2), params(f))
# end


text:
using CodecLz4
using Test

@testset "CodecLz4.jl" begin
    include("headers/lz4.jl")
    include("headers/lz4frame.jl")
    include("headers/lz4hc.jl")
    include("frame_compression.jl")
    include("hc_compression.jl")
    include("lz4_compression.jl")
    include("simple_compression.jl")
end


text:
@testset "Constants" begin
    @test Markers.SOI === 0xFFD8
end

@testset "Constructors" begin
    @test_throws ErrorException Marker(0xCAFE)
    @test Marker(0xFF01) === Marker{0xFF01}()
end

@testset "Operators" begin
    @test Marker{:SOF1}() !== Marker{0xFFC1}()
    @test Marker{:SOF1}() == Marker{0xFFC1}()
end

@testset "Traits" begin
    @test sizeof(Marker{:JPG0}) == 2
    @test sizeof(Marker{:JPG0}()) == 2
end

@testset "Display" begin
    @test shown_string(Marker{:SOI}) == "Marker{:SOI}"
    @test shown_string(LowLevelJPEG.SOI) == "Marker{:SOI}()"
end

@testset "Functions" begin
    @test written_hex(LowLevelJPEG.EOI) == "ffd9"

    @test reinterpret(UInt16, Marker{:SOF0}()) == 0xFFC0
end


text:
fdir = "J:\\ENG\\PAC_Group\\CurrentProjects\\Electric ServoMotor Brake\\Plan & Proof\\Plan\\2019-9-13 Coil Winding"
_savefig = (plt,fname)-> begin
    fdir = "J:\\ENG\\PAC_Group\\CurrentProjects\\Electric ServoMotor Brake\\Plan & Proof\\Plan\\2019-9-13 Coil Winding"
    savefig(plt,joinpath(fdir,fname*".pdf"))
    savefig(plt,joinpath(fdir,fname*".png"))
end

plt=let
    step=1u"mm"
    id=70u"mm":step:140u"mm"
    od=80u"mm":step:150u"mm"
    len=26u"mm"
    power=50.0u"W"
    global f = x -> begin
        coil=CoilGeometry(x.id,x.od,x.len)
        if CoilFunctions.isvalidcoil(coil)
            CoilFunctions.minimumturns(coil,power,1.0)
        else
            NaN
        end
    end
    @time y=tmap(f,((id=id,od=od,len=len) for id in id, od in od))
    contour(ustrip.(id),ustrip.(od),y';fill=true,
        title="SBE7 24V Minimum Equivalent Turns",xlabel="ID (mm)",ylabel="OD (mm)",dpi=300)
    plot!([0, 100,100],[140,140,0];color=:red,linestyle=:dash,
        ylims=ustrip.(extrema(od)),xlims=ustrip.(extrema(id)),
        legend=:bottomright,label="SBE7 Reference")
end
_savefig(plt,"SBE7 24V Minimum Equivalent Turns")

plt=let
    step=1u"mm"
    id=70u"mm":step:140u"mm"
    od=80u"mm":step:150u"mm"
    len=26u"mm"
    power=50.0u"W"
    global f = x -> begin
        coil=CoilGeometry(x.id,x.od,x.len)
        if CoilFunctions.isvalidcoil(coil)
            CoilFunctions.minimumturns(coil,power,1.0;op=extrema)
        else
            (NaN, NaN)
        end
    end
    @time y=tmap(f,((id=id,od=od,len=len) for id in id, od in od))
    contour(ustrip.(id),ustrip.(od),map(x->x[2]-x[1],permutedims(y,(2,1)));fill=true,
        title="SBE7 24V Equivalent Range of Turn Count",xlabel="ID (mm)",ylabel="OD (mm)",dpi=300)
    plot!([0, 100,100],[140,140,0];color=:red,linestyle=:dash,
        ylims=ustrip.(extrema(od)),xlims=ustrip.(extrema(id)),
        legend=:bottomright,label="SBE7 Reference")
end
_savefig(plt,"SBE7 24V Equivalent Range of Turn Count")

plt=let
    step=0.1u"mm"
    id=70u"mm":step:140u"mm"
    od=80u"mm":step:150u"mm"
    len=26u"mm"
    power=50.0u"W"
    v = 24.0u"V"
    Ω = v^2/power |> u"Ω"
    global f = x -> begin
        coil=CoilGeometry(x.id,x.od,x.len)
        if CoilFunctions.isvalidcoil(coil)
            optimalcoil(coil,Ω,1.0).turns
        else
            NaN
        end
    end
    @time y=tmap(f,((id=id,od=od,len=len) for id in id, od in od))
    contour(ustrip.(id),ustrip.(od),y';fill=true,
        title="SBE7 24V Coil Minimum Turns",xlabel="ID (mm)",ylabel="OD (mm)",dpi=300)
    plot!([0, 100,100],[140,140,0];color=:red,linestyle=:dash,
        ylims=ustrip.(extrema(od)),xlims=ustrip.(extrema(id)),
        legend=:bottomright,label="SBE7 Reference")
end
_savefig(plt,"SBE7 24V Coil Minimum Turns")

plt=let
    step=0.1u"mm"
    id=80u"mm":step:120u"mm"
    od=130u"mm":step:150u"mm"
    len=26u"mm"
    power=50.0u"W"
    v = 24.0u"V"
    Ω = v^2/power |> u"Ω"
    global f = x -> begin
        coil=CoilGeometry(x.id,x.od,x.len)
        if CoilFunctions.isvalidcoil(coil)
            optimalcoil(coil,Ω,1.0).gauge
        else
            NaN
        end
    end
    @time y=tmap(f,((id=id,od=od,len=len) for id in id, od in od))
    contour(ustrip.(id),ustrip.(od),y';fill=true,
        title="SBE7 24V Coil Wire Gauge",xlabel="ID (mm)",ylabel="OD (mm)",dpi=300,
        levels=length(unique(y)))
end
_savefig(plt,"SBE7 24V Coil Wire Gauge")

max_stars_repo_path: scripts/Graphs for A3 #1.jl

max_issues_repo_path: scripts/Graphs for A3 #1.jl

max_forks_repo_path: scripts/Graphs for A3 #1.jl


text:
using MetaGraphs
using LightGraphs
using CSV

enum(vec) = Dict(v => i for (i, v) ∈ enumerate(vec))

@testset "Lattice graph" begin
   m = 4
   n = 4
   t = 4
   L = 128

   instance = "$(@__DIR__)/instances/chimera_droplets/$(L)power/001.txt"

   ig = ising_graph(instance)

   fg = factor_graph(ig, 2, cluster_assignment_rule=super_square_lattice((m, n, 2*t))    )

   @test collect(vertices(fg)) == [(i, j) for i ∈ 1:m for j ∈ 1:n]

   clv = []
   cle = []
   rank = rank_vec(ig)

   for v ∈ vertices(fg)
      cl = get_prop(fg, v, :cluster)
      push!(clv, vertices(cl))
      push!(cle, collect(edges(cl)))

      @test rank_vec(cl) == get_prop.(Ref(ig), vertices(cl), :rank)
   end

   # Check if graph is factored correctly
   @test isempty(intersect(clv...))
   @test isempty(intersect(cle...))
end

@testset "Factor graph builds on pathological instance" begin
   m = 3
   n = 4
   t = 3
   L = n * m * t

   instance = "$(@__DIR__)/instances/pathological/test_$(m)_$(n)_$(t).txt"

   ising = CSV.File(instance, types=[Int, Int, Float64], header=0, comment = "#")

   couplings = Dict((i, j) => v for (i, j, v) ∈ ising)

   cedges = Dict(
      ((1, 1), (1, 2)) => [(1, 4), (1, 5), (1, 6)],
      ((1, 1), (2, 1)) => [(1, 13)],
      ((1, 2), (1, 3)) => [(4, 7), (5, 7), (6, 8), (6, 9)],
      ((1, 2), (2, 2)) => [(6, 16), (6, 18), (5, 16)],
      ((2, 1), (2, 2)) => [(13, 16), (13, 18)],
      ((2, 2), (3, 2)) => [(18, 28)],
      ((3, 2), (3, 3)) => [(28, 31), (28, 32), (28, 33), (29, 31), (29, 32), (29, 33), (30, 31), (30, 32), (30, 33)]
   )

   cells = Dict(
      (1, 1) => [1],
      (1, 2) => [4, 5, 6],
      (1, 3) => [7, 8, 9],
      (1, 4) => [],
      (2, 1) => [13],
      (2, 2) => [16, 18],
      (2, 3) => [],
      (2, 4) => [],
      (3, 1) => [],
      (3, 2) => [28, 29, 30],
      (3, 3) => [31, 32, 33],
      (3, 4) => []
   )

   d = 2
   rank = Dict(
      c => fill(d, length(idx))
      for (c,idx) ∈ cells if !isempty(idx)
   )

   bond_dimensions = [2, 2, 8, 4, 2, 2, 8]

   fg = factor_graph(
      ising_graph(instance),
      spectrum=full_spectrum,
      cluster_assignment_rule=super_square_lattice((m, n, t)),
   )

   for (bd, e) in zip(bond_dimensions, edges(fg))
      pl, en, pr = get_prop(fg, e, :pl), get_prop(fg, e, :en), get_prop(fg, e, :pr)
      @test minimum(size(en)) == bd
   end

   for ((i, j), cedge) ∈ cedges
      pl, en, pr = get_prop(fg, i, j, :pl), get_prop(fg, i, j, :en), get_prop(fg, i, j, :pr)
      base_i = all_states(rank[i])
      base_j = all_states(rank[j])

      idx_i = enum(cells[i])
      idx_j = enum(cells[j])

      # Change it to test if energy is calculated using passed 'energy' function
      energy = zeros(prod(rank[i]), prod(rank[j]))

      for (ii, σ) ∈ enumerate(base_i)
         for (jj, η) ∈ enumerate(base_j)
            eij = 0.
            for (k, l) ∈ values(cedge)
               kk, ll = enum(cells[i])[k], enum(cells[j])[l]
               s, r = σ[idx_i[k]], η[idx_j[l]]
               J = couplings[k, l]
               eij += s * J * r
            end
            energy[ii, jj] = eij
         end
      end
      @test energy ≈ pl * (en * pr)
   end

   @testset "each cluster comprises expected cells" begin
   for v ∈ vertices(fg)
      cl = get_prop(fg, v, :cluster)

      @test issetequal(vertices(cl), cells[v])
   end
   end

   @testset "each edge comprises expected bunch of edges from source Ising graph" begin
   for e ∈ edges(fg)
      outer_edges = get_prop(fg, e, :outer_edges)

      @test issetequal(cedges[(src(e), dst(e))], [(src(oe), dst(oe)) for oe ∈ outer_edges])
   end
   end
end

@testset "Rank reveal correctly decomposes energy row-wise" begin
   energy = [[1 2 3]; [0 -1 0]; [1 2 3]]
   P, E = rank_reveal(energy, :PE)
   @test size(P) == (3, 2)
   @test size(E) == (2, 3)
   @test P * E ≈ energy
end

@testset "Rank reveal correctly decomposes energy column-wise" begin
   energy = [[1, 2, 3] [0, -1, 1] [1, 2, 3]]
   E, P = rank_reveal(energy, :EP)
   @test size(P) == (2, 3)
   @test size(E) == (3, 2)
   @test E * P ≈ energy
end

@testset "Rank reveal correctly decomposes energy into projector, energy, projector" begin
   #energy = [[1 2 3]; [0 -1 0]; [1 2 3]]
   energy = [
      [1.0  -0.5   1.5   0.0   0.0  -1.5   0.5  -1.0];
      [0.0   0.5   0.5   1.0  -1.0  -0.5  -0.5   0.0];
      [0.5   0.0   1.0   0.5  -0.5  -1.0   0.0  -0.5];
      [-0.5   1.0   0.0   1.5  -1.5   0.0  -1.0   0.5];
      [0.5  -1.0   0.0  -1.5   1.5   0.0   1.0  -0.5];
      [-0.5   0.0  -1.0  -0.5   0.5   1.0   0.0   0.5];
      [0.0  -0.5  -0.5  -1.0   1.0   0.5   0.5   0.0];
      [-1.0   0.5  -1.5   0.0   0.0   1.5  -0.5   1.0]
   ]
   Pl, E_old = rank_reveal(energy, :PE)
   @test size(Pl) == (8, 8)
   @test size(E_old) == (8, 8)
   @test Pl * E_old ≈ energy

   E, Pr = rank_reveal(E_old, :EP)
   @test size(Pr) == (8, 8)
   @test size(E) == (8, 8)
   @test E * Pr ≈ E_old
   @test Pl * E * Pr ≈ energy
end

function create_example_factor_graph()
   J12 = -1.0
   h1 = 0.5
   h2 = 0.75

   D = Dict((1, 2) => J12, (1, 1) => h1, (2, 2) => h2)
   ig = ising_graph(D)

   factor_graph(
      ig,
      Dict((1, 1) => 2, (1, 2) => 2),
      spectrum = full_spectrum,
      cluster_assignment_rule = Dict(1 => (1, 1), 2 => (1, 2)),
  )
end

const fg_state_to_spin = [
   ([1, 1], [-1, -1]), ([1, 2], [-1, 1]), ([2, 1], [1, -1]), ([2, 2], [1, 1])
]

@testset "Decoding solution gives correct spin assignment" begin
   fg = create_example_factor_graph()
   for (state, spin_values) ∈ fg_state_to_spin
      d = decode_factor_graph_state(fg, state)
      states = collect(values(d))[collect(keys(d))]
      @test states == spin_values
   end
end

"""
Instance below looks like this:

1 -- 2 -- 3
|    |    |
4 -- 5 -- 6
|    |    |
7 -- 8 -- 9

And we group the following spins together: [1, 2, 4, 5], [3, 6], [7, 8], [9].
"""
function create_larger_example_factor_graph()
   instance = Dict(
      (1, 1) => 0.5,
      (2, 2) => 0.25,
      (3, 3) => 0.3,
      (4, 4) => 0.1,
      (5, 5) => 0.0,
      (6, 6) => -2.0,
      (7, 7) => -1.0,
      (8, 8) => 2.0,
      (9, 9) => 3.1,
      (1, 2) => -1.0,
      (2, 3) => 1.0,
      (4, 5) => 0.5,
      (5, 6) => -0.3,
      (7, 8) => 0.1,
      (8, 9) => 2.2,
      (1, 4) => -1.7,
      (4, 7) => 1.2,
      (2, 5) => 0.2,
      (5, 8) => 0.3,
      (3, 6) => 1.1,
      (6, 9) => 0.7
   )

   ig = ising_graph(instance)

   assignment_rule = Dict(
      1 => (1, 1),
      2 => (1, 1),
      4 => (1, 1),
      5 => (1, 1),
      3 => (1, 2),
      6 => (1, 2),
      7 => (2, 1),
      8 => (2, 1),
      9 => (2, 2)
   )

   fg = factor_graph(
      ig,
      Dict{NTuple{2, Int}, Int}(),
      spectrum = full_spectrum,
      cluster_assignment_rule = assignment_rule,
   )

   ig, fg
end

function factor_graph_energy(fg, state)
   # This is highly inefficient, but simple, which makes it suitable for testing.
   # If such a function is needed elsewhere, we need to implement it properly.
   total_en = 0.0

   # Collect local terms from each cluster
   for (s, v) ∈ zip(state, vertices(fg))
      total_en += get_prop(fg, v, :spectrum).energies[s]
   end

   # Collect inter-cluster terms
   for edge ∈ edges(fg)
      i, j = fg.reverse_label_map[src(edge)], fg.reverse_label_map[dst(edge)]
      pl, en, pr = get_prop(fg, edge, :pl), get_prop(fg, edge, :en), get_prop(fg, edge, :pr)
      edge_energy = pl * en * pr
      total_en += edge_energy[state[i], state[j]]
   end

   total_en
end

@testset "Decoding solution gives spins configuration with corresponding energies" begin
   ig, fg = create_larger_example_factor_graph()

   # Corresponding bases sizes for each cluster are 16, 4, 4, 2.
   all_states = [[i, j, k, l] for i ∈ 1:16 for j ∈ 1:4 for k ∈ 1:4 for l ∈ 1:2]

   for state ∈ all_states
      d = decode_factor_graph_state(fg, state)
      spins = zeros(length(d))
      for (k, v) ∈ d
         spins[k] = v
      end
      @test factor_graph_energy(fg, state) ≈ energy(spins, ig)
   end
end


text:
"""
Propagates the ComesAfter constraint:
	it removes the rule from the domain if the predecessors sequence is in the ancestors
"""
function propagate(c::ComesAfter, context::GrammarContext, domain::Vector{Int})
	ancestors = get_rulesequence(context.originalExpr, context.nodeLocation[begin:end-1])  # remove the current node from the node sequence
	if c.rule in domain  # if rule is in domain, check the ancestors
		if containedin(c.predecessors, ancestors)
			return domain
		else
			return filter(e -> e != c.rule, domain)
		end
	else # if it is not in the domain, just return domain
		return domain
	end
end

"""
	Propagates Ordered constraint:
		removes every element from domain that does not have a necessary predecessor in the left subtree
"""
function propagate(c::Ordered, context::GrammarContext, domain::Vector{Int})
	expr = context.originalExpr
	rules_on_left = rulesonleft(context.originalExpr, context.nodeLocation)

	last_rule_index = 0
	for r in c.order
		r in rules_on_left ? last_rule_index = r : break
	end

	rules_to_remove = Set(c.order[last_rule_index+2:end]) # +2 because the one after the last index can be used

	return filter((x) -> !(x in rules_to_remove), domain)
end


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule Graphene_jll
using Base
using Base: UUID
import JLLWrappers

JLLWrappers.@generate_main_file_header("Graphene")
JLLWrappers.@generate_main_file("Graphene", UUID("75302f13-0b7e-5bab-a6d1-23fa92e4c2ea"))
end  # module Graphene_jll


text:
#http://en.wikipedia.org/wiki/Stochastic_matrix

I = zeros(4,4);
[I[i,i] = 1 for i in 1:4];

f = open("m3.txt","r")
T = readdlm(f,',');
close(f);

Ep = [0 1 0 0]*inv(I - T)*[1,1,1,1];

println("Expected lifetime for the mouse is $Ep")


text:
module TestReadDatastores

using Test, FASTX, ReadDatastores, BioSequences

include("paired-reads.jl")
include("long-reads.jl")
include("linked-reads.jl")

@testset "Error Messages" begin
    buf = IOBuffer()

    Base.showerror(buf, ReadDatastores.MissingMagicError("myreads.prseq"))
    @test String(take!(buf)) == "MissingMagicError: the file myreads.prseq does not appear to be a valid read datastore file, it does not begin with the expected magic bytes."

    Base.showerror(buf, ReadDatastores.DatastoreTypeError{PairedReads{DNAAlphabet{2}}}("myreads.prseq", ReadDatastores.LongDS))
    @test String(take!(buf)) == "DatastoreTypeError: myreads.prseq contains a long read datastore and cannot be opened as a ReadDatastores.PairedReads{BioSequences.DNAAlphabet{2}}"

    Base.showerror(buf, ReadDatastores.DatastoreVersionError{PairedReads{DNAAlphabet{2}}}(UInt16(2)))
    @test String(take!(buf)) == "DatastoreVersionError: file format version of paired read datastore file (v2) is deprecated: this version of ReadDatastores.jl supports v$(Int(ReadDatastores.PairedDS_Version))"

    Base.showerror(buf, ReadDatastores.DatastoreVersionError{LongReads{DNAAlphabet{2}}}(UInt16(2)))
    @test String(take!(buf)) == "DatastoreVersionError: file format version of long read datastore file (v2) is deprecated: this version of ReadDatastores.jl supports v$(Int(ReadDatastores.LongDS_Version))"

    Base.showerror(buf, ReadDatastores.DatastoreVersionError{LinkedReads{DNAAlphabet{2}}}(UInt16(2)))
    @test String(take!(buf)) == "DatastoreVersionError: file format version of linked read datastore file (v2) is deprecated: this version of ReadDatastores.jl supports v$(Int(ReadDatastores.LinkedDS_Version))"

    Base.showerror(buf, ReadDatastores.DatastoreEncodingError{PairedReads{DNAAlphabet{2}}}("myreads.prseq", 4))
    @test String(take!(buf)) == "DatastoreEncodingError: myreads.prseq encodes reads using 4 bits per element and cannot be opened as a ReadDatastores.PairedReads{BioSequences.DNAAlphabet{2}}"
end

end # module


text:
using DynamicalSystemsBase

function eom_ar1_unidir(x, p, n)
    a₁, b₁, c_xy, σ = (p...,)
    x, y = (x...,)
    ξ₁ = rand(Normal(0, σ))
    ξ₂ = rand(Normal(0, σ))

    dx = a₁*x + ξ₁
    dy = b₁*y + c_xy*x + ξ₂
    return SVector{2}(dx, dy)
end

function ar1_unidir(;uᵢ = rand(2), a₁ = 0.90693, b₁ = 0.40693, c_xy = 0.5, σ = 0.40662)
    p = [a₁, b₁, c_xy, σ]
    DiscreteDynamicalSystem(eom_ar1_unidir, uᵢ, p)
end

vars = (1, 2)
npts, tstep = 50, 50
d_xind = Uniform(2.5, 5.5)
d_yind = Uniform(2.5, 5.5)
d_xval = Uniform(0.01, 0.2)
d_yval = Uniform(0.01, 0.2)

X, Y = example_uncertain_indexvalue_datasets(ar1_unidir(c_xy = 0.5), npts, vars, tstep = tstep,
    d_xind = d_xind, d_yind = d_yind,
    d_xval = d_xval, d_yval = d_yval);

time_grid = -20:100:2540
n_draws = 10000 # draws per uncertain value
n_bins = length(time_grid) - 1
wts = rand(length(X))

# Values in each bin represented as RawValues
b = BinnedWeightedResampling(RawValues, time_grid, wts, n_draws)
bc, vs = bin(Y, b);

@test vs isa Vector{Vector{T}} where T
@test length(vs) == n_bins

# Values in each bin represented as UncertainScalarKDE
b_kde = BinnedWeightedResampling(UncertainScalarKDE, time_grid, wts, n_draws)
Y_binned = bin(Y, b_kde);

@test Y_binned isa AbstractUncertainIndexValueDataset

# Values in each bin represented as UncertainScalarPopulation
b_pop = BinnedWeightedResampling(UncertainScalarPopulation, time_grid, wts, n_draws)
Y_binned = bin(Y, b_pop);

@test Y_binned isa AbstractUncertainIndexValueDataset


text:
using QuantumGate
using LinearAlgebra

# 初期状態
# Notation: |0> = [1, 0], |1> = [0, 1]
#
firstbit = [1; 0]
ψ = [1; 0]
ket = firstbit ⊗ ψ

# アダマールテストでHadを測定する回路
res = Hadamardgate(1, 2)*controlled_Ugate(1, 2, 2, Had)*Hadamardgate(1, 2)*ket

# 1つ目のqubitが0である確率は
# p = (1 + Re<ψ|H|ψ>)/2
# if |ψ> = |0>, p = (1 + 1/√2)/2
measured = (P⁰ ⊗ σ⁰) * res

dot(measured, measured)
(1 + 1/sqrt(2))/2

# 測定後の状態はHの固有状態である。
Had * measured[1:2] ≈ measured[1:2]

# 1つ目のqubitが1である確率は
# p = (1 - Re<ψ|H|ψ>)/2
# if |ψ> = |0>, p = (1 - 1/√2)/2
measured = (P¹ ⊗ σ⁰) * res

dot(measured, measured)
(1 - 1/sqrt(2))/2

# 測定後の状態はHの固有状態である。
Had * measured[3:4] ≈ -measured[3:4]


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

mutable struct PublicIPAddressSku2 <: SwaggerModel
    name::Any # spec type: Union{ Nothing, String } # spec name: name

    function PublicIPAddressSku2(;name=nothing)
        o = new()
        validate_property(PublicIPAddressSku2, Symbol("name"), name)
        setfield!(o, Symbol("name"), name)
        o
    end
end # type PublicIPAddressSku2

const _property_map_PublicIPAddressSku2 = Dict{Symbol,Symbol}(Symbol("name")=>Symbol("name"))
const _property_types_PublicIPAddressSku2 = Dict{Symbol,String}(Symbol("name")=>"String")
Base.propertynames(::Type{ PublicIPAddressSku2 }) = collect(keys(_property_map_PublicIPAddressSku2))
Swagger.property_type(::Type{ PublicIPAddressSku2 }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_PublicIPAddressSku2[name]))}
Swagger.field_name(::Type{ PublicIPAddressSku2 }, property_name::Symbol) =  _property_map_PublicIPAddressSku2[property_name]

const _allowed_PublicIPAddressSku2_name = ["Basic", "Standard"]

function check_required(o::PublicIPAddressSku2)
    true
end

function validate_property(::Type{ PublicIPAddressSku2 }, name::Symbol, val)
    if name === Symbol("name")
        Swagger.validate_param(name, "PublicIPAddressSku2", :enum, val, _allowed_PublicIPAddressSku2_name)
    end
end


text:
#Subject
"""
    setkelrange!(data::T, range::ElimRange{:point}; kelauto = false) where T <: PKSubject

Set `range` for subject `data`. Set `kelauto` if possible.
"""
function setkelrange!(data::T, range::ElimRange{:point}; kelauto = false) where T <: PKSubject
    if range.kelend > length(data) throw(ArgumentError("Kel endpoint out of range")) end
    data.kelrange = range
    setkelauto!(data, kelauto)
    data
end
#DS ind Int
"""
    setdosetime!(data::DataSet{T}, dosetime::DoseTime, ind::Int) where T <: PKSubject
"""
function setkelrange!(data::DataSet{T}, range::ElimRange{:point}, ind::Int; kelauto = false) where T <: PKSubject
    setkelrange!(data[ind], range; kelauto = kelauto)
    data
end
#DS iter Int
"""
    setkelrange!(data::DataSet{T}, range::ElimRange{:point}, inds::Union{Vector{Int}, UnitRange{Int}, Tuple{Vararg{Int}}}; kelauto = false)
"""
function setkelrange!(data::DataSet{T}, range::ElimRange{:point}, inds::Union{Vector{Int}, UnitRange{Int}, Tuple{Vararg{Int}}}; kelauto = false) where T <: PKSubject
    for i in inds
        setkelrange!(data[i], range; kelauto = kelauto)
    end
    data
end
#DS all
"""
    setkelrange!(data::DataSet{T}, range::ElimRange{:point}; kelauto = false) where T <: PKSubject
"""
function setkelrange!(data::DataSet{T}, range::ElimRange{:point}; kelauto = false) where T <: PKSubject
    for i = 1:length(data)
        setkelrange!(data[i], range; kelauto = kelauto)
    end
    data
end
#DS Dict
"""
    setkelrange!(data::DataSet{T}, range::ElimRange{:point}, sort::Dict; kelauto = false) where T <: PKSubject
"""
function  setkelrange!(data::DataSet{T}, range::ElimRange{:point}, sort::Dict; kelauto = false) where T <: PKSubject
    for i = 1:length(data)
        if sort ⊆ data[i].id setkelrange!(data[i], range; kelauto = kelauto) end
    end
    data
end
#GET subj
"""
    getkelrange(data::T) where T <: PKSubject
"""
function getkelrange(data::T) where T <: PKSubject
    data.kelrange
end


text:
export make_lattice, makelattice
export Lattice
export dimension

"""
    Lattice{S, O}

Represent a lattice.

# Arguments
* `unitcell::UnitCell{S, O}`
* `hypercube::Hypercube`
* `bravais_coordinates::Vector{Vector{Int}}`
* `supercell::UnitCell{Tuple{S, Vector{Int}}, Tuple{O, Vector{Int}}}`
"""
struct Lattice{S, O}
    unitcell::UnitCell{S, O}
    hypercube::Hypercube
    bravais_coordinates::Vector{Vector{Int}}
    supercell::UnitCell{Tuple{S, Vector{Int}}, Tuple{O, Vector{Int}}}
end

function Base.:(==)(lhs::Lattice{S, O}, rhs::Lattice{S, O}) where {S, O}
    return (
        lhs.unitcell == rhs.unitcell
        && lhs.hypercube == rhs.hypercube
        && lhs.bravais_coordinates == rhs.bravais_coordinates
    )
end

"""
    makelattice(unitcell, shape)

Create a `Lattice` from `unitcell` and `shape`.
Generating translations are fouond using `find_generators`.

# Arguments
- `unitcell`
- `shape`: shape of the cluster
"""
function makelattice(
    unitcell::UnitCell{S, O},
    shape_matrix::AbstractMatrix{<:Integer}
) where {S, O}
    dim = dimension(unitcell)
    if size(shape_matrix) != (dim, dim)
        throw(DimensionMismatch("unitcell and shape_matrix should have the same dimension"))
    end

    hypercube = Hypercube(shape_matrix)
    generator_translations = find_generators(hypercube)
    unitcell_coordinates = generate_coordinates(hypercube, generator_translations)

    new_latticevectors = unitcell.latticevectors * hypercube.shape_matrix

    new_unitcell = makeunitcell(new_latticevectors; SiteType=Tuple{S, Vector{Int}}, OrbitalType=Tuple{O, Vector{Int}})
    for uc in unitcell_coordinates
        for (sitename, sitecoord) in unitcell.sites
            cc = fract2carte(unitcell, sitecoord)
            new_cc = cc + unitcell.latticevectors * uc
            new_sitecoord = carte2fract(new_unitcell, new_cc)
            new_sitename = (sitename, uc)
            addsite!(new_unitcell, new_sitename, new_sitecoord)
        end
        for (orbitalname, orbitalsiteindex) in unitcell.orbitals
            sitename = getsitename(unitcell, orbitalsiteindex)
            new_orbitalname = (orbitalname, uc)
            new_sitename = (sitename, uc)
            new_orbitalsiteindex = getsiteindex(new_unitcell, new_sitename)
            addorbital!(new_unitcell, new_orbitalname, new_orbitalsiteindex)
        end
    end
    return Lattice{S, O}(unitcell, hypercube, unitcell_coordinates, new_unitcell)
end
#TODO unit testing for lattice with wrong dimensions

"""
    makelattice(unitcell, shape, [generators])

Create a lattice with periodic boundary condition, using the unitcell, shape, and
translation generators.

# Arguments
* `unitcell::UnitCell{O}`
* `shape_matrix::AbstractMatrix{<:Integer}`: shape of the Bravais lattice. This can also be
a single integer, which is equivalent to a identity matrix times the number.

# Optional Arguments
* `generator_translations::AbstractMatrix{<:Integer}`
"""
function makelattice(
    unitcell::UnitCell{S, O},
    shape_matrix::AbstractMatrix{<:Integer},
    generator_translations::AbstractMatrix{<:Integer},
) where {S, O}
    dim = dimension(unitcell)
    if size(shape_matrix) != (dim, dim)
        throw(DimensionMismatch("unitcell and shape_matrix should have the same dimension"))
    end
    hypercube = Hypercube(shape_matrix)
    unitcell_coordinates = generate_coordinates(hypercube, generator_translations)

    new_latticevectors = unitcell.latticevectors * hypercube.shape_matrix

    new_unitcell = makeunitcell(new_latticevectors; SiteType=Tuple{S, Vector{Int}}, OrbitalType=Tuple{O, Vector{Int}})
    for uc in unitcell_coordinates
        for (sitename, sitecoord) in unitcell.sites
            cc = fract2carte(unitcell, sitecoord)
            new_cc = cc + unitcell.latticevectors * uc
            new_sitecoord = carte2fract(new_unitcell, new_cc)
            new_sitename = (sitename, uc)
            addsite!(new_unitcell, new_sitename, new_sitecoord)
        end
        for (orbitalname, orbitalsiteindex) in unitcell.orbitals
            sitename = getsitename(unitcell, orbitalsiteindex)
            new_orbitalname = (orbitalname, uc)
            new_sitename = (sitename, uc)
            new_orbitalsiteindex = getsiteindex(new_unitcell, new_sitename)
            addorbital!(new_unitcell, new_orbitalname, new_orbitalsiteindex)
        end
    end

    return Lattice{S, O}(unitcell, hypercube, unitcell_coordinates, new_unitcell)
end

"""
    makelattice(unitcell, scale)

Construct a `Lattice` with `unitcell` having shape `scale`.
Works for one-dimensional lattice.
"""
function makelattice(unitcell::UnitCell, scale::Integer)
    return makelattice(unitcell, hcat(scale), hcat(1))
end

const make_lattice = makelattice

"""
    dimension(lattice)

Spatial dimension of the lattice.

# Arguments
- `lattice::Lattice`
"""
dimension(lattice::Lattice) = dimension(lattice.unitcell)


text:
using FITSIO, Plots

f_strip = FITS("weather_STRIP_ERA5.fits")
f_alma  = FITS("weather_file/weather_ALMA.fits")
f_qubic = FITS("weather_file/weather_QUBIC_3.dits")
f_reijo = FITS("weather_file/weather_Atacama.fits")

seasons_matrix_strip = Array{Float64, 2}(undef, 12, 24)
seasons_matrix_qubic = Array{Float64, 2}(undef, 12, 24)
seasons_matrix_reijo = Array{Float64, 2}(undef, 12, 24)

for i in 2:13
    for j in 1:24
        v = read(f_strip[i], "TQV")[50, j]
        v_q = read(f_qubic[i], "TQV")[50, j]
        v_r = read(f_reijo[i], "TQV")[50, j]
        seasons_matrix_strip[i-1,j] = v #* (1000/997)
        seasons_matrix_qubic[i-1,j] = v_q #* (1000/997)
        seasons_matrix_reijo[i-1,j] = v_r
    end
end

gradient=ColorGradient([:black, :aqua])

heatmap(seasons_matrix_strip, xlabel="Hour", ylabel="Month", zlabel="TQV [mm]", ytick=1:12, xticks=1:24, title = "TQV fluctuations CL 95%", c=gradient, clims=(0,20))

heatmap(seasons_matrix_reijo, xlabel="Hour", ylabel="Month", zlabel="TQV [mm]", ytick=1:12, xticks=1:24, title = "TQV fluctuations CL 50%", c=gradient, clims=(0,9))

# heatmap(read(f_ten[10], "TQV"))

max_stars_repo_name: algebrato/atmosphere

max_issues_repo_name: algebrato/atmosphere

max_forks_repo_name: algebrato/atmosphere


text:
# UnitCommitment.jl: Optimization Package for Security-Constrained Unit Commitment
# Copyright (C) 2020, UChicago Argonne, LLC. All rights reserved.
# Released under the modified BSD license. See COPYING.md for more details.

using UnitCommitment
using JuMP
import UnitCommitment:
    ArrCon2000,
    CarArr2006,
    DamKucRajAta2016,
    Formulation,
    Gar1962,
    KnuOstWat2018,
    MorLatRam2013,
    PanGua2016,
    XavQiuWanThi2019

if ENABLE_LARGE_TESTS
    using Gurobi
end

function _small_test(formulation::Formulation)::Nothing
    instances = ["matpower/case118/2017-02-01", "test/case14"]
    for instance in instances
        # Should not crash
        UnitCommitment.build_model(
            instance = UnitCommitment.read_benchmark(instance),
            formulation = formulation,
        )
    end
    return
end

function _large_test(formulation::Formulation)::Nothing
    instances = ["pglib-uc/ca/Scenario400_reserves_1"]
    for instance in instances
        instance = UnitCommitment.read_benchmark(instance)
        model = UnitCommitment.build_model(
            instance = instance,
            formulation = formulation,
            optimizer = Gurobi.Optimizer,
        )
        UnitCommitment.optimize!(
            model,
            XavQiuWanThi2019.Method(two_phase_gap = false, gap_limit = 0.1),
        )
        solution = UnitCommitment.solution(model)
        @test UnitCommitment.validate(instance, solution)
    end
    return
end

function _test(formulation::Formulation)::Nothing
    _small_test(formulation)
    if ENABLE_LARGE_TESTS
        _large_test(formulation)
    end
end

@testset "formulations" begin
    _test(Formulation())
    _test(Formulation(ramping = ArrCon2000.Ramping()))
    # _test(Formulation(ramping = DamKucRajAta2016.Ramping()))
    _test(
        Formulation(
            ramping = MorLatRam2013.Ramping(),
            startup_costs = MorLatRam2013.StartupCosts(),
        ),
    )
    _test(Formulation(ramping = PanGua2016.Ramping()))
    _test(Formulation(pwl_costs = Gar1962.PwlCosts()))
    _test(Formulation(pwl_costs = CarArr2006.PwlCosts()))
    _test(Formulation(pwl_costs = KnuOstWat2018.PwlCosts()))
end


text:
include("common.jl")

example = """
5483143223
2745854711
5264556173
6141336146
6357385478
4167524645
2176841721
6882881134
4846848554
5283751526
"""

function to_matrix(input::AbstractString)
    rows = Vector{Int64}[]
    for line in split(input, '\n')
        isempty(line) && continue
        push!(rows, [parse(Int64, x) for x in line])
    end
    return transpose(reduce(hcat, rows)) |> collect
end

function apply_flashes!(x::AbstractMatrix{Int64}, flashed::AbstractMatrix{Bool})
    @assert size(x) == size(flashed)
    ni = size(flashed, 1)
    nj = size(flashed, 2)
    @inbounds for i in 1:ni
        for j in 1:nj
            flashed[i, j] || continue

            # Update neighbours.
            for oi in max(1, (i - 1)):min((i + 1), ni)
                for oj in max(1, (j - 1)):min((j + 1), nj)
                    # Skip the thing that has flashed.
                    oi == i && oj == j && continue
                    x[oi, oj] += 1
                end
            end
        end
    end
    return x
end

function step!(x::AbstractMatrix{Int64})
    x .+= 1  # Energy level increases.

    flashed = zeros(Bool, size(x))
    while true
        # An octopus shouldn't flash more than once.
        this_flashed = @. ((x > 9) & (!flashed))
        sum(this_flashed) == 0 && break  # Done flashing
        apply_flashes!(x, this_flashed)
        flashed .|= this_flashed  # Update the record of all flashes.
    end

    x[flashed] .= 0
    return sum(flashed)
end

function evolve!(x::AbstractMatrix, num_steps)
    total = 0
    for _ in 1:num_steps
        total += step!(x)
    end
    return total
end

function get_ans(input::AbstractString, n=100)
    x = to_matrix(input)
    return evolve!(x, n)
end

@assert get_ans(example, 10) == 204
@assert get_ans(example) == 1656
get_ans(load_input(11))

# Part 2

function get_ans2(input::AbstractString)
    x = to_matrix(input)
    step = 1
    while true
        num_flashes = step!(x)
        num_flashes == length(x) && return step
        step += 1
    end
end

@assert get_ans2(example) == 195
get_ans2(load_input(11))


text:
module semiclassical

using QuantumOpticsBase
import Base: ==
import ..timeevolution: integrate, recast!, jump, integrate_mcwf, jump_callback, as_vector, QO_CHECKS
import LinearAlgebra: normalize, normalize!

using Random, LinearAlgebra
import OrdinaryDiffEq

# TODO: Remove imports
import DiffEqCallbacks, RecursiveArrayTools.copyat_or_push!
Base.@pure pure_inference(fout,T) = Core.Compiler.return_type(fout, T)

using ..timeevolution

const QuantumState{B} = Union{Ket{B}, DenseOperator{B,B}}
const DecayRates = Union{Nothing, Vector{Float64}, Matrix{Float64}}

"""
Semi-classical state.

It consists of a quantum part, which is either a `Ket` or a `DenseOperator` and
a classical part that is specified as a complex vector of arbitrary length.
"""
mutable struct State{B<:Basis,T<:QuantumState{B},C<:Vector{ComplexF64}}
    quantum::T
    classical::C
    function State(quantum::T, classical::C) where {B<:Basis,T<:QuantumState{B},C<:Vector{ComplexF64}}
        new{B,T,C}(quantum, classical)
    end
end

Base.length(state::State) = length(state.quantum) + length(state.classical)
Base.copy(state::State) = State(copy(state.quantum), copy(state.classical))
normalize!(state::State{B,T}) where {B,T<:Ket} = normalize!(state.quantum)
normalize(state::T) where {B,K<:Ket,T<:State{B,K}} = State(normalize(state.quantum),copy(state.classical))

function ==(a::State, b::State)
    QuantumOpticsBase.samebases(a.quantum, b.quantum) &&
    length(a.classical)==length(b.classical) &&
    (a.classical==b.classical) &&
    (a.quantum==b.quantum)
end

QuantumOpticsBase.expect(op, state::State) = expect(op, state.quantum)
QuantumOpticsBase.variance(op, state::State) = variance(op, state.quantum)
QuantumOpticsBase.ptrace(state::State, indices::Vector{Int}) = State(ptrace(state.quantum, indices), state.classical)

QuantumOpticsBase.dm(x::State{B,T}) where {B<:Basis,T<:Ket{B}} = State(dm(x.quantum), x.classical)

"""
    semiclassical.schroedinger_dynamic(tspan, state0, fquantum, fclassical[; fout, ...])

Integrate time-dependent Schrödinger equation coupled to a classical system.

# Arguments
* `tspan`: Vector specifying the points of time for which the output should
        be displayed.
* `psi0`: Initial semi-classical state [`semiclassical.State`](@ref).
* `fquantum`: Function `f(t, psi, u) -> H` returning the time and or state
        dependent Hamiltonian.
* `fclassical`: Function `f(t, psi, u, du)` calculating the possibly time and
        state dependent derivative of the classical equations and storing it
        in the vector `du`.
* `fout=nothing`: If given, this function `fout(t, state)` is called every time
        an output should be displayed. ATTENTION: The given state is neither
        normalized nor permanent!
* `kwargs...`: Further arguments are passed on to the ode solver.
"""
function schroedinger_dynamic(tspan, state0::S, fquantum::Function, fclassical::Function;
                fout::Union{Function,Nothing}=nothing,
                kwargs...) where {B<:Basis,T<:Ket{B},S<:State{B,T}}
    tspan_ = convert(Vector{Float64}, tspan)
    dschroedinger_(t::Float64, state::S, dstate::S) = dschroedinger_dynamic(t, state, fquantum, fclassical, dstate)
    x0 = Vector{ComplexF64}(undef, length(state0))
    recast!(state0, x0)
    state = copy(state0)
    dstate = copy(state0)
    integrate(tspan_, dschroedinger_, x0, state, dstate, fout; kwargs...)
end

"""
    semiclassical.master_dynamic(tspan, state0, fquantum, fclassical; <keyword arguments>)

Integrate time-dependent master equation coupled to a classical system.

# Arguments
* `tspan`: Vector specifying the points of time for which output should
        be displayed.
* `rho0`: Initial semi-classical state [`semiclassical.State`](@ref).
* `fquantum`: Function `f(t, rho, u) -> (H, J, Jdagger)` returning the time
        and/or state dependent Hamiltonian and Jump operators.
* `fclassical`: Function `f(t, rho, u, du)` calculating the possibly time and
        state dependent derivative of the classical equations and storing it
        in the complex vector `du`.
* `fout=nothing`: If given, this function `fout(t, state)` is called every time
        an output should be displayed. ATTENTION: The given state is not
        permanent!
* `kwargs...`: Further arguments are passed on to the ode solver.
"""
function master_dynamic(tspan, state0::State{B,T}, fquantum, fclassical;
                rates::DecayRates=nothing,
                fout::Union{Function,Nothing}=nothing,
                tmp::T=copy(state0.quantum),
                kwargs...) where {B<:Basis,T<:DenseOperator{B,B}}
    tspan_ = convert(Vector{Float64}, tspan)
    function dmaster_(t::Float64, state::S, dstate::S) where {B<:Basis,T<:DenseOperator{B,B},S<:State{B,T}}
        dmaster_h_dynamic(t, state, fquantum, fclassical, rates, dstate, tmp)
    end
    x0 = Vector{ComplexF64}(undef, length(state0))
    recast!(state0, x0)
    state = copy(state0)
    dstate = copy(state0)
    integrate(tspan_, dmaster_, x0, state, dstate, fout; kwargs...)
end

function master_dynamic(tspan, state0::State{B,T}, fquantum, fclassical; kwargs...) where {B<:Basis,T<:Ket{B}}
    master_dynamic(tspan, dm(state0), fquantum, fclassical; kwargs...)
end

"""
    semiclassical.mcwf_dynamic(tspan, psi0, fquantum, fclassical, fjump_classical; <keyword arguments>)

Calculate MCWF trajectories coupled to a classical system.

# Arguments
* `tspan`: Vector specifying the points of time for which output should
        be displayed.
* `psi0`: Initial semi-classical state [`semiclassical.State`](@ref) featuring
        a `Ket`(@ref).
* `fquantum`: Function `f(t, psi, u) -> (H, J, Jdagger)` returning the time
        and/or state dependent Hamiltonian and Jump operators.
* `fclassical`: Function `f(t, psi, u, du)` calculating the possibly time and
        state dependent derivative of the classical equations and storing it
        in the complex vector `du`.
* `fjump_classical`: Function `f(t, psi, u, i)` performing a classical jump when a
        quantum jump of the i-th jump operator occurs.
* `fout=nothing`: If given, this function `fout(t, state)` is called every time
        an output should be displayed. ATTENTION: The given state is not
        permanent!
* `display_beforeevent`: Choose whether or not an additional point should be saved
        before a jump occurs. Default is false.
* `display_afterevent`: Choose whether or not an additional point should be saved
        after a jump occurs. Default is false.
* `display_jumps=false`: If set to true, an additional list of times and indices
        is returned. These correspond to the times at which a jump occured and
        the index of the jump operators with which the jump occured, respectively.
* `kwargs...`: Further arguments are passed on to the ode solver.
"""
function mcwf_dynamic(tspan, psi0::State{B,T}, fquantum, fclassical, fjump_classical;
                seed=rand(UInt),
                rates::DecayRates=nothing,
                fout::Union{Function,Nothing}=nothing,
                kwargs...) where {B<:Basis,T<:Ket{B}}
    tspan_ = convert(Vector{Float64}, tspan)
    tmp=copy(psi0.quantum)
    function dmcwf_(t::Float64, psi::S, dpsi::S) where {B<:Basis,T<:Ket{B},S<:State{B,T}}
        dmcwf_h_dynamic(t, psi, fquantum, fclassical, rates, dpsi, tmp)
    end
    j_(rng, t::Float64, psi, psi_new) = jump_dynamic(rng, t, psi, fquantum, fclassical, fjump_classical, psi_new, rates)
    x0 = Vector{ComplexF64}(undef, length(psi0))
    recast!(psi0, x0)
    psi = copy(psi0)
    dpsi = copy(psi0)
    integrate_mcwf(dmcwf_, j_, tspan_, psi, seed, fout; kwargs...)
end

function recast!(state::State{B,T,C}, x::C) where {B<:Basis,T<:QuantumState{B},C<:Vector{ComplexF64}}
    N = length(state.quantum)
    copyto!(x, 1, state.quantum.data, 1, N)
    copyto!(x, N+1, state.classical, 1, length(state.classical))
    x
end

function recast!(x::C, state::State{B,T,C}) where {B<:Basis,T<:QuantumState{B},C<:Vector{ComplexF64}}
    N = length(state.quantum)
    copyto!(state.quantum.data, 1, x, 1, N)
    copyto!(state.classical, 1, x, N+1, length(state.classical))
end

function dschroedinger_dynamic(t::Float64, state::State{B,T}, fquantum::Function,
            fclassical::Function, dstate::State{B,T}) where {B<:Basis,T<:Ket{B}}
    fquantum_(t, psi) = fquantum(t, state.quantum, state.classical)
    timeevolution.dschroedinger_dynamic(t, state.quantum, fquantum_, dstate.quantum)
    fclassical(t, state.quantum, state.classical, dstate.classical)
end

function dmaster_h_dynamic(t::Float64, state::State{B,T}, fquantum::Function,
            fclassical::Function, rates::DecayRates, dstate::State{B,T}, tmp::T) where {B<:Basis,T<:DenseOperator{B,B}}
    fquantum_(t, rho) = fquantum(t, state.quantum, state.classical)
    timeevolution.dmaster_h_dynamic(t, state.quantum, fquantum_, rates, dstate.quantum, tmp)
    fclassical(t, state.quantum, state.classical, dstate.classical)
end

function dmcwf_h_dynamic(t::Float64, psi::T, fquantum::Function, fclassical::Function, rates::DecayRates,
                    dpsi::T, tmp::K) where {T,K}
    fquantum_(t, rho) = fquantum(t, psi.quantum, psi.classical)
    timeevolution.dmcwf_h_dynamic(t, psi.quantum, fquantum_, rates, dpsi.quantum, tmp)
    fclassical(t, psi.quantum, psi.classical, dpsi.classical)
end

function jump_dynamic(rng, t::Float64, psi::T, fquantum::Function, fclassical::Function, fjump_classical::Function, psi_new::T, rates::DecayRates) where T<:State
    result = fquantum(t, psi.quantum, psi.classical)
    QO_CHECKS[] && @assert 3 <= length(result) <= 4
    J = result[2]
    if length(result) == 3
        rates_ = rates
    else
        rates_ = result[4]
    end
    i = jump(rng, t, psi.quantum, J, psi_new.quantum, rates_)
    fjump_classical(t, psi_new.quantum, psi.classical, i)
    psi_new.classical .= psi.classical
    return i
end

function jump_callback(jumpfun::Function, seed, scb, save_before!::Function,
                        save_after!::Function, save_t_index::Function, psi0::State)
    tmp = copy(psi0)
    psi_tmp = copy(psi0)
    rng = MersenneTwister(convert(UInt, seed))
    jumpnorm = Ref(rand(rng))
    n = length(psi0.quantum)
    djumpnorm(x::Vector{ComplexF64}, t::Float64, integrator) = norm(x[1:n])^2 - (1-jumpnorm[])

    function dojump(integrator)
        x = integrator.u
        t = integrator.t

        affect! = scb.affect!
        save_before!(affect!,integrator)
        recast!(x, psi_tmp)
        i = jumpfun(rng, t, psi_tmp, tmp)
        recast!(tmp, x)
        save_after!(affect!,integrator)
        save_t_index(t,i)

        jumpnorm[] = rand(rng)
        return nothing
    end

    return OrdinaryDiffEq.ContinuousCallback(djumpnorm,dojump,
                     save_positions = (false,false))
end
as_vector(psi::State{B,K}) where {B,K<:Ket} = [psi.quantum.data; psi.classical]

end # module


text:
module Wire

using ..Ahorn, Maple

const placements = Ahorn.PlacementDict(
    "Wire" => Ahorn.EntityPlacement(
        Maple.Wire,
        "line"
    )
)

Ahorn.nodeLimits(entity::Maple.Wire) = 1, 1

function Ahorn.selection(entity::Maple.Wire)
    nodes = get(entity.data, "nodes", ())
    x, y = Ahorn.position(entity)

    res = Ahorn.Rectangle[Ahorn.Rectangle(x - 4, y - 4, 8, 8)]

    for node in nodes
        nx, ny = node

        push!(res, Ahorn.Rectangle(nx - 4, ny - 4, 8, 8))
    end

    return res
end

wireColorHex = "595866"
wireColor = (89, 88, 102, 1) ./ (255, 255, 255, 1)
wireColorSelected = (Ahorn.colors.selection_selected_fc)

function renderWire(ctx::Ahorn.Cairo.CairoContext, entity::Maple.Wire, color::Ahorn.colorTupleType=wireColor)
    x, y = Ahorn.position(entity)

    start = (x, y)
    stop = get(entity.data, "nodes", [start])[1]
    control = (start .+ stop) ./ 2 .+ (0, 24)

    curve = Ahorn.SimpleCurve(start, stop, control)
    Ahorn.drawSimpleCurve(ctx, curve, color, thickness=1)
end

Ahorn.renderSelectedAbs(ctx::Ahorn.Cairo.CairoContext, entity::Maple.Wire) = renderWire(ctx, entity, wireColorSelected)

function Ahorn.renderAbs(ctx::Ahorn.Cairo.CairoContext, entity::Maple.Wire, room::Maple.Room)
    # Make sure Alpha is 1
    rawColor = Ahorn.argb32ToRGBATuple(parse(Int, get(entity.data, "color", wireColorHex), base=16))[1:3] ./ 255
    color = (rawColor..., 1.0)

    renderWire(ctx, entity, color)
end

end


text:
# Testing the 4-10 particles example interval [0,0.6]
# Static Parameters
f_c = 20
x_max = 1.0
num_x = 10*f_c
filter_x = ones(2*f_c+1)

# Dynamic Parameters
K  = 2
tau = 1.0/(K*1.0)
v_max = 0.5
num_v = 10

# Variability of weights
minWeights = 0.9
maxWeights = 1.1

# Variability of number of particles
min_number_part = 4
max_number_part = 10

# Number of generated examples
num_trials = 1000

# Iterations of minimization algorithm
algIter = 100;

# test case
test_case = () -> TestCases.cloud_1d_full(x_max, v_max, minWeights, maxWeights, K, tau, rand(min_number_part:max_number_part))
# WARNING: changing anything in the test_case will generate troubles if the rejection_sampling algorithm
# in Utils.js is not either updated for it, or disabled.

# Noise parameters
noises_data = linspace(0,0.1,5)
noises_data = noises_data[2:end]
noises_position = linspace(0,0.06,5)
noises_position = noises_position[2:end]

### Location of data folder
dataFolder = "data/1Dsimulations"
# loaded data for rejection sampling of simulations
bins = np.load(dataFolder*"/separationDistribBins_06interval.npy")
density = np.load(dataFolder*"/separationDistribVal_06interval.npy")

println(" WARNING: The rejection sampling algorithm is being used, with interval of interest: [",bins[1]," ",bins[end],"]. To change it check the separationDistribution.jl file.")

### Cases to be tested

noiseless_dynamic = true
noiseless_static = true
noise_dynamic = true
noise_static = true
curvature_static = true

cases = (noiseless_dynamic, noiseless_static, noise_dynamic, noise_static, curvature_static)


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

mutable struct IoK8sApiAuthorizationV1SubjectRulesReviewStatus <: SwaggerModel
    evaluationError::Any # spec type: Union{ Nothing, String } # spec name: evaluationError
    incomplete::Any # spec type: Union{ Nothing, Bool } # spec name: incomplete
    nonResourceRules::Any # spec type: Union{ Nothing, Vector{IoK8sApiAuthorizationV1NonResourceRule} } # spec name: nonResourceRules
    resourceRules::Any # spec type: Union{ Nothing, Vector{IoK8sApiAuthorizationV1ResourceRule} } # spec name: resourceRules

    function IoK8sApiAuthorizationV1SubjectRulesReviewStatus(;evaluationError=nothing, incomplete=nothing, nonResourceRules=nothing, resourceRules=nothing)
        o = new()
        validate_property(IoK8sApiAuthorizationV1SubjectRulesReviewStatus, Symbol("evaluationError"), evaluationError)
        setfield!(o, Symbol("evaluationError"), evaluationError)
        validate_property(IoK8sApiAuthorizationV1SubjectRulesReviewStatus, Symbol("incomplete"), incomplete)
        setfield!(o, Symbol("incomplete"), incomplete)
        validate_property(IoK8sApiAuthorizationV1SubjectRulesReviewStatus, Symbol("nonResourceRules"), nonResourceRules)
        setfield!(o, Symbol("nonResourceRules"), nonResourceRules)
        validate_property(IoK8sApiAuthorizationV1SubjectRulesReviewStatus, Symbol("resourceRules"), resourceRules)
        setfield!(o, Symbol("resourceRules"), resourceRules)
        o
    end
end # type IoK8sApiAuthorizationV1SubjectRulesReviewStatus

const _property_map_IoK8sApiAuthorizationV1SubjectRulesReviewStatus = Dict{Symbol,Symbol}(Symbol("evaluationError")=>Symbol("evaluationError"), Symbol("incomplete")=>Symbol("incomplete"), Symbol("nonResourceRules")=>Symbol("nonResourceRules"), Symbol("resourceRules")=>Symbol("resourceRules"))
const _property_types_IoK8sApiAuthorizationV1SubjectRulesReviewStatus = Dict{Symbol,String}(Symbol("evaluationError")=>"String", Symbol("incomplete")=>"Bool", Symbol("nonResourceRules")=>"Vector{IoK8sApiAuthorizationV1NonResourceRule}", Symbol("resourceRules")=>"Vector{IoK8sApiAuthorizationV1ResourceRule}")
Base.propertynames(::Type{ IoK8sApiAuthorizationV1SubjectRulesReviewStatus }) = collect(keys(_property_map_IoK8sApiAuthorizationV1SubjectRulesReviewStatus))
Swagger.property_type(::Type{ IoK8sApiAuthorizationV1SubjectRulesReviewStatus }, name::Symbol) = Union{Nothing,eval(Base.Meta.parse(_property_types_IoK8sApiAuthorizationV1SubjectRulesReviewStatus[name]))}
Swagger.field_name(::Type{ IoK8sApiAuthorizationV1SubjectRulesReviewStatus }, property_name::Symbol) =  _property_map_IoK8sApiAuthorizationV1SubjectRulesReviewStatus[property_name]

function check_required(o::IoK8sApiAuthorizationV1SubjectRulesReviewStatus)
    (getproperty(o, Symbol("incomplete")) === nothing) && (return false)
    (getproperty(o, Symbol("nonResourceRules")) === nothing) && (return false)
    (getproperty(o, Symbol("resourceRules")) === nothing) && (return false)
    true
end

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


text:
function ParallelPark(method=:none)

    opts = SolverOptions(
        cost_tolerance_intermediate=1e-3,
        active_set_tolerance=1e-4
    )

    # model
    model = RobotZoo.DubinsCar()
    n,m = size(model)
    N = 101
    tf = 3.

    # cost
    x0 = @SVector [0., 0., 0.]
    xf = @SVector [0., 1., 0.]
    Qf = 100.0*Diagonal(@SVector ones(n))
    Q = (1e-2)*Diagonal(@SVector ones(n))
    R = (1e-2)*Diagonal(@SVector ones(m))

    # constraints
    u_bnd = 2.
    x_min = @SVector [-0.25, -0.001, -Inf]
    x_max = @SVector [0.25, 1.001, Inf]
    bnd = BoundConstraint(n,m,x_min=x_min,x_max=x_max,u_min=-u_bnd,u_max=u_bnd)
    goal = GoalConstraint(xf)

    # Constraint vals
    con_bnd = ConstraintVals(bnd, 1:N-1)
    con_goal = ConstraintVals(goal, N:N)

    # problem
    U = [@SVector fill(0.1,m) for k = 1:N-1]
    obj = LQRObjective(Q,R,Qf,xf,N)

    conSet = ConstraintSet(n,m,[con_bnd, con_goal], N)

    prob = Problem(model, obj, xf, tf, constraints=conSet, x0=x0, U0=U)

    return prob, opts
end


text:
# Decoding is done by using Viterbi Algorithm
# Computes in polynomial time

"""
Scores for the first tag in the tagging sequence.
"""
function preds_first(c::CRF, y)
    c.W[c.n + 1, onecold(y, 1:length(y))]
end

"""
Scores for the last tag in the tagging sequence.
"""
function preds_last(c::CRF, y)
    c.W[onecold(y, 1:length(y)), c.n + 2]
end

"""
Scores for the tags other than the starting one.
"""
function preds_single(c::CRF, y, y_prev)
    c.W[onecold(y_prev, 1:length(y_prev)), onecold(y, 1:length(y))]
end

# Helper for forward pass, returns max_probs and corresponding arg_max for all the labels
function forward_pass_unit(k)
    α_idx = [i[1]  for i in argmax(k, dims=1)]
    α = [k[j, i] for (i,j) in enumerate(α_idx)]
    return α, α_idx
end

"""
Computes the forward pass for viterbi algorithm.
"""
function _decode(c::CRF, x, init_vit_vars)
    α_idx = zeros(Int, c.n + 2, length(x))

    forward_var, α_idx[:, 1] = forward_pass_unit((c.W .+ x[1]') .+ init_vit_vars)

    for i in 2:length(x)
        forward_var, α_idx[:, i] = forward_pass_unit((c.W .+ x[i]') .+ forward_var')
    end

    labels = zeros(Int, length(x))
    labels[end] = argmax(forward_var + (c.W[:, c.n + 2])')[2]

    for i in reverse(2:length(x))
        labels[i - 1] =  α_idx[labels[i], i]
    end

    @assert α_idx[labels[1], 1] == c.n + 1 # Check for START Tag
    return onehotseq(labels, c.n)
end

onehotseq(seq, num_labels) = [onehot(i, 1:num_labels) for i in seq]

"""
    viterbi_decode(::CRF, input_sequence)

Predicts the most probable label sequence of `input_sequence`.
"""
function viterbi_decode(c::CRF, x_seq, init_vit_vars)
    length(x_seq) == 0 && throw("Input sequence is empty")
    return _decode(cpu(c), cpu.(x_seq), cpu(init_vit_vars))
end

# function predict(c::CRF, x_seq)
#     viterbi_decode(c, x_seq)
# end


text:
@doc raw"""
Models.CorrelatedGaussian(ndims)

Creates a highly-correlated Gaussian with the given dimensionality.

```math
\mathbf\theta \sim \mathcal{N}\left(2\mathbf{1}, \mathbf{I}\right)
```
```math
\Sigma_{ij} = \begin{cases} 1 &\quad i=j \\ 0.95 &\quad i\neq j \end{cases}
```
```math
\mathcal{L}(\mathbf\theta) = \mathcal{N}\left(\mathbf\theta | \mathbf{0}, \mathbf\Sigma \right)
```

the analytical evidence of the model is

```math
Z = \mathcal{N}\left(2\mathbf{1} | \mathbf{0}, \mathbf\Sigma + \mathbf{I} \right)
```

## Examples
```jldoctest
julia> model, lnZ = Models.CorrelatedGaussian(10);

julia> lnZ
-12.482738597926607
```
"""
function CorrelatedGaussian(ndims)
    priors = fill(Normal(2, 1), ndims)
    Σ = fill(0.95, ndims, ndims)
    Σ[diagind(Σ)] .= 1
    cent_dist = MvNormal(Σ)
    loglike(X) = logpdf(cent_dist, X)

    model = NestedModel(loglike, priors)
    true_lnZ = logpdf(MvNormal(fill(2, ndims), Σ + I), zeros(ndims))
    return model, true_lnZ
end


text:
require( bytestring(ENV["OCCA_DIR"], "/lib/occa.jl") )

entries = 5

Serial_Info   = "mode = Serial"
OpenMP_Info   = "mode = OpenMP  , schedule = compact, chunk = 10";
OpenCL_Info   = "mode = OpenCL  , platformID = 0, deviceID = 0";
CUDA_Info     = "mode = CUDA    , deviceID = 0";
Pthreads_Info = "mode = Pthreads, threadCount = 4, schedule = compact, pinnedCores = [0, 0, 1, 1]";
COI_Info      = "mode = COI     , deviceID = 0";

device = occa.device(Serial_Info);

a  = Float32[1 - i for i in 1:entries]
b  = Float32[i     for i in 1:entries]
ab = Float32[0     for i in 1:entries]

o_a  = occa.malloc(device, a);
o_b  = occa.malloc(device, b);
o_ab = occa.malloc(device, ab);

addVectors = occa.buildKernelFromSource(device,
                                        "addVectors.okl",
                                        "addVectors")

occa.runKernel(addVectors,
               (entries, Int32),
               o_a, o_b, o_ab)

occa.memcpy(ab, o_ab)

println(ab)


text:
module SmartPlayTest

using  Test

import Cribbage
import Cribbage.SmartPlay: SmartPlayer
import Cribbage.TwoPlayer
import Cribbage.ThreePlayer
import Cribbage.FourPlayer

@testset "Smart Players" begin
    p₁ = SmartPlayer("one")
    @test p₁.name == "one"
    p₂ = SmartPlayer("two")
    @test p₂.name == "two"
    @test p₁ ≠ p₂
end

@testset "Smart Gameplay" begin

    # Logging.with_logger(logger) do
        p₁ = SmartPlayer("player one")
        p₂ = SmartPlayer("player two")
        game = Cribbage.TwoPlayer.TwoPlayerGame(p₁, p₂)
        gamestate = Cribbage.TwoPlayer.TwoPlayerGameState(
            game, rand([p₁, p₂])
        )

        gamestates = Cribbage.playgame(gamestate)

        p₃ = SmartPlayer("player three")
        game = Cribbage.ThreePlayer.ThreePlayerGame(p₁, p₂, p₃)
        gamestate = Cribbage.ThreePlayer.ThreePlayerGameState(
            game, rand([p₁, p₂, p₃])
        )
        gamestates = Cribbage.playgame(gamestate)

        p₄ = SmartPlayer("player four")
        game = Cribbage.FourPlayer.FourPlayerGame(p₁, p₂, p₃, p₄)
        gamestate = Cribbage.FourPlayer.FourPlayerGameState(
            game, rand([p₁, p₂, p₃, p₄])
        )
        gamestates = Cribbage.playgame(gamestate)
    # end

end

end  # module SmartPlayTest


text:
class C {
  int{L} x = 7;
  void m() {
      int{H} y = 42;
      this.x = ({L}) y;
  }
    public static void main(String[] args) {
	new C().m();
    }
}

max_stars_repo_name: HarvardPL/cryptoerase

max_issues_repo_name: HarvardPL/cryptoerase

max_forks_repo_name: HarvardPL/cryptoerase


text:
using TerminalSpinners
using Test

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


text:
function eratosphen(n)
    ser=fill(true,n)
    ser[1]=false
    k=2
    while k !== nothing && k<n
        ser[k^2:k:end] .= false
        k=findnext(ser, k+1)
    end
    return findall(ser)
end


text:
module Bridges

using Compat

using MathOptInterface
const MOI = MathOptInterface
const MOIU = MOI.Utilities

const MOIB = MOI.Bridges # used in macro

const SVF = MOI.SingleVariable
const VVF = MOI.VectorOfVariables
const SAF{T} = MOI.ScalarAffineFunction{T}
const VAF{T} = MOI.VectorAffineFunction{T}
const SQF{T} = MOI.ScalarQuadraticFunction{T}
const VQF{T} = MOI.VectorQuadraticFunction{T}

const VI = MOI.VariableIndex
const CI = MOI.ConstraintIndex

include("bridge.jl")
include("bridgeoptimizer.jl")
include("singlebridgeoptimizer.jl")
include("lazybridgeoptimizer.jl")

# This is used by JuMP and removes the need to update JuMP everytime a bridge is added
MOIU.@model(AllBridgedConstraints,
            (),
            (MOI.EqualTo, MOI.LessThan, MOI.GreaterThan, MOI.Interval,),
            (MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone,
             MOI.RotatedSecondOrderCone, MOI.GeometricMeanCone,
             MOI.PositiveSemidefiniteConeSquare,
             MOI.LogDetConeTriangle, MOI.RootDetConeTriangle),
            (),
            (MOI.SingleVariable,),
            (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction),
            (MOI.VectorOfVariables,),
            (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction))
"""
    full_bridge_optimizer(model::MOI.ModelLike, ::Type{T}) where T

Returns a `LazyBridgeOptimizer` bridging `model` for every bridge defined in
this package and for the coefficient type `T`.
"""
function full_bridge_optimizer(model::MOI.ModelLike, ::Type{T}) where T
    cache = MOIU.UniversalFallback(AllBridgedConstraints{T}())
    bridged_model = LazyBridgeOptimizer(model, cache)
    add_bridge(bridged_model, GreaterToLessBridge{T})
    add_bridge(bridged_model, LessToGreaterBridge{T})
    add_bridge(bridged_model, NonnegToNonposBridge{T})
    add_bridge(bridged_model, NonposToNonnegBridge{T})
    add_bridge(bridged_model, ScalarizeBridge{T})
    add_bridge(bridged_model, VectorizeBridge{T})
    add_bridge(bridged_model, ScalarSlackBridge{T})
    add_bridge(bridged_model, VectorSlackBridge{T})
    add_bridge(bridged_model, SplitIntervalBridge{T})
    add_bridge(bridged_model, QuadtoSOCBridge{T})
    add_bridge(bridged_model, GeoMeanBridge{T})
    add_bridge(bridged_model, SquarePSDBridge{T})
    add_bridge(bridged_model, LogDetBridge{T})
    add_bridge(bridged_model, RootDetBridge{T})
    add_bridge(bridged_model, RSOCBridge{T})
    add_bridge(bridged_model, RSOCtoPSDBridge{T})
    add_bridge(bridged_model, SOCtoPSDBridge{T})
    return bridged_model
end

include("flip_sign_bridge.jl")
@bridge GreaterToLess GreaterToLessBridge () (MOI.GreaterThan,) () () (MOI.SingleVariable,) (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction) () ()
@bridge LessToGreater LessToGreaterBridge () (MOI.LessThan,) () () (MOI.SingleVariable,) (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction) () ()
@bridge NonnegToNonpos NonnegToNonposBridge () () (MOI.Nonnegatives,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction)
@bridge NonposToNonneg NonposToNonnegBridge () () (MOI.Nonpositives,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction)
include("vectorizebridge.jl")
@bridge Vectorize VectorizeBridge () (MOI.EqualTo, MOI.LessThan, MOI.GreaterThan,) () () (MOI.SingleVariable,) (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction) () ()
include("scalarizebridge.jl")
@bridge Scalarize ScalarizeBridge () () (MOI.Zeros, MOI.Nonnegatives, MOI.Nonpositives) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction)
include("slackbridge.jl")
@bridge ScalarSlack ScalarSlackBridge () (MOI.Interval, MOI.LessThan, MOI.GreaterThan) () () () (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction) () ()
@bridge(VectorSlack, VectorSlackBridge,  (), (),
    (MOI.Nonnegatives, MOI.Nonpositives, MOI.SecondOrderCone,
    MOI.RotatedSecondOrderCone, MOI.GeometricMeanCone,
    MOI.PositiveSemidefiniteConeSquare, MOI.PositiveSemidefiniteConeTriangle, MOI.LogDetConeTriangle,
    MOI.RootDetConeTriangle),
    (MOI.PowerCone, MOI.DualPowerCone, MOI.SOS1, MOI.SOS2), (), (), (),
    (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction)
    )
include("intervalbridge.jl")
@bridge SplitInterval SplitIntervalBridge () (MOI.Interval,) () () (MOI.SingleVariable,) (MOI.ScalarAffineFunction, MOI.ScalarQuadraticFunction) () ()
include("rsocbridge.jl")
@bridge RSOC RSOCBridge () () (MOI.RotatedSecondOrderCone,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)
include("quadtosocbridge.jl")
@bridge QuadtoSOC QuadtoSOCBridge () (MOI.LessThan, MOI.GreaterThan) () () () (MOI.ScalarQuadraticFunction,) () ()
include("geomeanbridge.jl")
@bridge GeoMean GeoMeanBridge () () (MOI.GeometricMeanCone,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)
include("squarepsdbridge.jl")
@bridge SquarePSD SquarePSDBridge () () (MOI.PositiveSemidefiniteConeSquare,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction, MOI.VectorQuadraticFunction)
include("detbridge.jl")
@bridge LogDet LogDetBridge () () (MOI.LogDetConeTriangle,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)
@bridge RootDet RootDetBridge () () (MOI.RootDetConeTriangle,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)
include("soctopsdbridge.jl")
@bridge SOCtoPSD SOCtoPSDBridge () () (MOI.SecondOrderCone,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)
@bridge RSOCtoPSD RSOCtoPSDBridge () () (MOI.RotatedSecondOrderCone,) () () () (MOI.VectorOfVariables,) (MOI.VectorAffineFunction,)

end # module


text:
"""
    count(::Type{MangalNetwork}, d::MangalDataset, query::Pair...)

Return the number of networks that belong to a dataset, according to an optional
series of queries.
"""
function count(::Type{MangalNetwork}, d::MangalDataset, query::Pair...)
    return count(MangalNetwork, "dataset_id" => d.id, query...)
end

"""
    count(::Type{MangalInteraction}, n::MangalNetwork, query::Pair...)

Return the number of interactions that belong to a network, according to an optional
series of queries.
"""
function count(::Type{MangalInteraction}, n::MangalNetwork, query::Pair...)
    return count(MangalInteraction, "network_id" => n.id, query...)
end

"""
    count(::Type{MangalNode}, n::MangalNetwork, query::Pair...)

Return the number of nodes that belong to a network, according to an optional
series of queries.
"""
function count(::Type{MangalNode}, n::MangalNetwork, query::Pair...)
    return count(MangalNode, "network_id" => n.id, query...)
end

"""
    count(::Type{MangalNode}, n::MangalReferenceTaxon, query::Pair...)

Return the number of nodes that are instances of a reference taxon, according to an optional
series of queries.
"""
function count(::Type{MangalNode}, r::MangalReferenceTaxon, query::Pair...)
    return count(MangalNode, "taxonomy_id" => r.id, query...)
end


text:
using ExplainableAI
using Documenter
using Literate

EXAMPLE_DIR = joinpath(@__DIR__, "literate")
OUT_DIR = joinpath(@__DIR__, "src/generated")

# Use Literate.jl to generate docs and notebooks of examples
for example in readdir(EXAMPLE_DIR)
    EXAMPLE = joinpath(EXAMPLE_DIR, example)

    Literate.markdown(EXAMPLE, OUT_DIR; documenter=true) # markdown for Documenter.jl
    Literate.notebook(EXAMPLE, OUT_DIR) # .ipynb notebook
    Literate.script(EXAMPLE, OUT_DIR) # .jl script
end

makedocs(;
    modules=[ExplainableAI],
    authors="Adrian Hill",
    sitename="ExplainableAI.jl",
    format=Documenter.HTML(; prettyurls=get(ENV, "CI", "false") == "true", assets=String[]),
    pages=[
        "Home" => "index.md",
        "Getting started" => "generated/example.md",
        "Advanced LRP" => "generated/advanced_lrp.md",
        "API Reference" => "api.md",
    ],
)

deploydocs(; repo="github.com/adrhill/ExplainableAI.jl")


text:
module MDToolbox

using FileIO
using NetCDF
using NLsolve
using FFTW
using Bio3DView
using CuArrays

# satndard library
using Printf
using Statistics
using LinearAlgebra
using Dates
using Base.Threads

# package code goes here
export KB_kcalpermol, KB_kjpermol
export AbstractTrajectory, TrjArray, select_atom
export readdcd, readnetcdf, writenetcdf, readpsf, writepsf, readpdb, writepdb, readcrd
export centerofmass, decenter, decenter!, superimpose, superimpose_serial, getrmsd, meanstructure, getrmsf, rotate, rotate!
export getdistance, getangle, getdihedral
export getmsd
export propagate_mcmc, propagate_md
export ksdensity, ksdensity_serial, getpmf
export wham, wham_iteration
export mbar
export sp_delta_pmf, sp_design_matrix, sp_design_matrix_atom, sp_lsquares, sp_admm, sp_descent, sp_standardize!, sp_standardize, sp_cumulate_pmf, sp_cumulate_pmf_atom
export afmize, AfmizeConfig, translateafm, getafmposterior, getafmposteriors_alpha
export Asd, readasd
export viewstruc, gpu, gpu32, cpu, cpu64

# constants
const KB_kcalpermol = 0.0019872041 #Boltzmann constant taken from wikipedia
const KB_kjpermol = 0.0083144621 #Boltzmann constant taken from wikipedia

# codes
include("trjarray.jl")
include("fileio.jl")
include("structure.jl")
include("dynamics.jl")
include("ksdensity.jl")
include("wham.jl")
include("mbar.jl")
include("sparse.jl")
include("afm.jl")
include("asd.jl")
include("utils.jl")

end # module


text:
"""
  iteratively_hypergraph_embedding(h::Hypergraph; dims=2, epoch=10000)

Embeds hypernodes and hyperedges into a `dims`-dimensional vector.
Time complexity is O(`epoch` * `dims` * ∑ᵥdeg(v)).

**Arguments**
* h : Hypergraph.
* dims : Dimension of the embedded vector.
* max_epoch : Number of iterations to be performed to obtain each vector.
As large as possible is better.
* ϵ : Convergence condition.
"""

function iteratively_hypergraph_embedding(h::Hypergraph, dims::Int=2,
    max_epoch::Int=10000, ϵ::Float64=1e-8, he_vec=nothing, hn_vec=nothing)
  @assert length(get_connected_components(h)) == 1 "Not Connected."
  N = nhv(h)
  M = nhe(h)
  if isnothing(he_vec) he_vec=rand(Uniform(-1, 1), dims, M) end
  if isnothing(hn_vec) hn_vec=rand(Uniform(-1, 1), dims, N) end
  foreach(normalize!, eachcol(he_vec))
  foreach(normalize!, eachcol(hn_vec))
  old_sum_dot = Inf

  for e in 1:max_epoch
    for d in 1:dims hn_vec[d, :] .-= mean(hn_vec[d, :]) end
    he_vec_new = zeros(dims, M)
    for he in 1:M
      hns = keys(getvertices(h, he))
      gᵣ = length(hns)
      for hn in hns
        dᵥ = length(gethyperedges(h, hn))
        he_vec_new[:, he] += hn_vec[:, hn]
      end
    end

    for d in 1:dims he_vec_new[d, :] .-= mean(he_vec_new[d, :]) end
    new_sum_dot = sum(dot.(eachcol(he_vec), eachcol(he_vec_new)))
    foreach(normalize!, eachcol(he_vec_new))
    he_vec = deepcopy(he_vec_new)

    for d in 1:dims he_vec[d, :] .-= mean(he_vec[d, :]) end
    hn_vec_new = zeros(dims, N)
    for hn in 1:N
      hes = keys(gethyperedges(h, hn))
      dᵥ = length(hes)
      for he in hes
        gᵣ = length(getvertices(h, he))
        hn_vec_new[:, hn] += he_vec[:, he]
      end
    end
    for d in 1:dims hn_vec_new[d, :] .-= mean(hn_vec_new[d, :]) end
    foreach(normalize!, eachcol(hn_vec_new))
    hn_vec = deepcopy(hn_vec_new)

    if abs(new_sum_dot - old_sum_dot) < ϵ
      break
    end
    old_sum_dot = new_sum_dot
  end

  return hn_vec, he_vec
end


text:
struct Error <: Exception
    str::String
    T
    code
end

Error(buf::AbstractString, T, code, pos, totallen) = Error(buf, T, code)
Error(buf::AbstractVector{UInt8}, T, code, pos, totallen) = Error(String(copy(buf)), T, code)

function Error(buf::IO, T, code, pos, totallen)
    fastseek!(buf, pos - 1)
    bytes = read(buf, totallen)
    return Error(String(bytes), T, code)
end

function Base.showerror(io::IO, e::Error)
    c = e.code
    println(io, "Parsers.Error ($(codes(c))):")
    println(io, text(c))
    println(io, "attempted to parse $(e.T) from: \"$(escape_string(e.str))\"")
end

"""
A bitmask value, with various bits corresponding to different parsing signals and scenarios.

`Parsers.xparse` returns a `code` value with various bits set according to the various scenarios
encountered while parsing a value.

    * `INVALID`: there are a number of invalid parsing states, all include the INVALID bit set (check via `Parsers.invalid(code)`)
    * `OK`: signals specifically that a valid value of type `T` was parsed (check via `Parsers.ok(code)`)
    * `SENTINEL`: signals that a valid sentinel value was detected while parsing, passed via the `sentinel` keyword argument to `Parsers.Options` (check via `Parsers.sentinel(code)`)
    * `QUOTED`: a `openquotechar` from `Parsers.Options` was detected at the beginning of parsing (check via `Parsers.quoted(code)`)
    * `DELIMITED`: a `delim` character or string from `Parsers.Options` was detected while parsing (check via `Parsers.delimited(code)`)
    * `NEWLINE`: a non-quoted newline character (`'\n'`), return character (`'\r'`), or CRLF (`"\r\n"`) was detected while parsing (check via `Parsers.newline(code)`)
    * `EOF`: the end of file was reached while parsing
    * `ESCAPED_STRING`: an `escapechar` from `Parsers.Options` was encountered while parsing (check via `Parsers.escapedstring(code)`)
    * `INVALID_QUOTED_FIELD`: a `openquotechar` were detected when parsing began, but no corresponding `closequotechar` were found to correctly close a quoted field, this is usually a fatal parsing error because parsing will continue until EOF to look for the close quote character (check via `Parsers.invalidquotedfield(code)`)
    * `INVALID_DELIMITER`: a `delim` character or string were eventually detected, but not at the expected position (directly after parsing a valid value), indicating there are extra, invalid characters between a valid value and the expected delimiter (check via `Parsers.invaliddelimiter(code)`)
    * `OVERFLOW`: overflow occurred while parsing an `Integer` value type (check via `Parsers.overflow(code)`)

One additional convenience function is provided, `Parsers.quotednotescaped(code)`, which checks if a value was quoted,
but didn't contain any escape characters, useful to indicate if a string may be used "as-is", instead of needing to be unescaped.
"""
const ReturnCode = Int16

const SUCCESS = 0b0000000000000000 % ReturnCode
const INVALID = 0b1000000000000000 % ReturnCode

# success flags
const OK                   = 0b0000000000000001 % ReturnCode
const SENTINEL             = 0b0000000000000010 % ReturnCode

# property flags
const QUOTED               = 0b0000000000000100 % ReturnCode
const DELIMITED            = 0b0000000000001000 % ReturnCode
const NEWLINE              = 0b0000000000010000 % ReturnCode
const EOF                  = 0b0000000000100000 % ReturnCode
const ESCAPED_STRING       = 0b0000001000000000 % ReturnCode

# invalid flags
const INVALID_QUOTED_FIELD = 0b1000000001000000 % ReturnCode
const INVALID_DELIMITER    = 0b1000000010000000 % ReturnCode
const OVERFLOW             = 0b1000000100000000 % ReturnCode

ok(x::ReturnCode) = (x & (OK | INVALID)) == OK
invalid(x::ReturnCode) = x < 0
sentinel(x::ReturnCode) = (x & SENTINEL) == SENTINEL
quoted(x::ReturnCode) = (x & QUOTED) == QUOTED
delimited(x::ReturnCode) = (x & DELIMITED) == DELIMITED
newline(x::ReturnCode) = (x & NEWLINE) == NEWLINE
escapedstring(x::ReturnCode) = (x & ESCAPED_STRING) == ESCAPED_STRING
invalidquotedfield(x::ReturnCode) = (x & INVALID_QUOTED_FIELD) == INVALID_QUOTED_FIELD
invaliddelimiter(x::ReturnCode) = (x & INVALID_DELIMITER) == INVALID_DELIMITER
overflow(x::ReturnCode) = (x & OVERFLOW) == OVERFLOW
quotednotescaped(x::ReturnCode) = (x & (QUOTED | ESCAPED_STRING)) == QUOTED

memcmp(a::Ptr{UInt8}, b::Ptr{UInt8}, len::Int) = ccall(:memcmp, Cint, (Ptr{UInt8}, Ptr{UInt8}, Csize_t), a, b, len) == 0

@inline function checksentinel(source::AbstractVector{UInt8}, pos, len, sentinel, debug)
    sentinelpos = 0
    startptr = pointer(source, pos)
    # sentinel is an iterable of Tuple{Ptr{UInt8}, Int}, sorted from longest sentinel string to shortest
    for (ptr, ptrlen) in sentinel
        if pos + ptrlen - 1 <= len
            match = memcmp(startptr, ptr, ptrlen)
            if match
                if debug
                    println("matched sentinel value: \"$(escape_string(unsafe_string(ptr, ptrlen)))\"")
                end
                sentinelpos = pos + ptrlen
                break
            end
        end
    end
    return sentinelpos
end

@inline function checksentinel(source::IO, pos, len, sentinel, debug)
    sentinelpos = 0
    origpos = position(source)
    for (ptr, ptrlen) in sentinel
        matched = match!(source, ptr, ptrlen)
        if matched
            if debug
                println("matched sentinel value: \"$(escape_string(unsafe_string(ptr, ptrlen)))\"")
            end
            sentinelpos = pos + ptrlen
            break
        end
    end
    fastseek!(source, origpos)
    return sentinelpos
end

@inline function checkdelim(source::AbstractVector{UInt8}, pos, len, (ptr, ptrlen))
    startptr = pointer(source, pos)
    delimpos = pos
    if pos + ptrlen - 1 <= len
        match = memcmp(startptr, ptr, ptrlen)
        if match
            delimpos = pos + ptrlen
        end
    end
    return delimpos
end

@inline function checkdelim(source::IO, pos, len, (ptr, ptrlen))
    delimpos = pos
    matched = match!(source, ptr, ptrlen)
    if matched
        delimpos = pos + ptrlen
    end
    return delimpos
end

@inline function match!(source::IO, ptr, ptrlen)
    b = peekbyte(source)
    c = unsafe_load(ptr)
    if b != c
        return false
    elseif ptrlen == 1
        incr!(source)
        return true
    end
    pos = position(source)
    i = 1
    while true
        incr!(source)
        if Base.eof(source)
            fastseek!(source, pos)
            return false
        end
        b = peekbyte(source)
        c = unsafe_load(ptr + i)
        if b != c
            fastseek!(source, pos)
            return false
        elseif i == ptrlen - 1
            break
        end
        i += 1
    end
    incr!(source)
    return true
end

"""
    Parsers.fastseek!(io::IO, n::Integer)

    Without valididty checks, seek an `IO` to desired byte position `n`. Used in Parsers.jl to
    seek back to a previous location already parsed.
"""
function fastseek! end

fastseek!(io::IO, n::Integer) = seek(io, n)
function fastseek!(io::IOBuffer, n::Integer)
    io.ptr = n + 1
    return
end
fastseek!(io::AbstractVector{UInt8}, n::Integer) = nothing

"""
    Parsers.readbyte(io::IO)::UInt8

    Consume a single byte from an `IO` without checking `eof(io)`.
"""
function readbyte end

"""
    Parsers.peekbyte(io::IO)::UInt8

    Return, but do not consume, the next byte from an `IO` without checking `eof(io)`.
"""
function peekbyte end

incr!(io::IO) = readbyte(io)
readbyte(from::IO) = Base.read(from, UInt8)
peekbyte(from::IO) = UInt8(Base.peek(from))

function readbyte(from::IOBuffer)
    i = from.ptr
    @inbounds byte = from.data[i]
    from.ptr = i + 1
    return byte
end

function peekbyte(from::IOBuffer)
    @inbounds byte = from.data[from.ptr]
    return byte
end

function incr!(from::IOBuffer)
    from.ptr += 1
    return
end

incr!(from::AbstractVector{UInt8}) = nothing
peekbyte(from::IO, pos) = peekbyte(from)
function peekbyte(from::AbstractVector{UInt8}, pos)
    @inbounds b = from[pos]
    return b
end

eof(source::AbstractVector{UInt8}, pos::Integer, len::Integer) = pos > len
eof(source::IO, pos::Integer, len::Integer) = Base.eof(source)

function text(r)
    str = ""
    if r & QUOTED > 0
        str = "encountered an opening quote character, initial value parsing "
    else
        str = "initial value parsing "
    end
    if r & OK > 0
        str *= "succeeded"
    else
        str *= "failed"
    end
    if r & (~INVALID & OVERFLOW) > 0
        str *= ", value overflowed"
    end
    if r & SENTINEL > 0
        str *= ", a sentinel value was parsed"
    end
    if r & ESCAPED_STRING > 0
        str *= ", encountered escape character"
    end
    if r & (~INVALID & INVALID_QUOTED_FIELD) > 0
        str *= ", invalid quoted value"
    end
    if r & DELIMITED > 0
        str *= ", a valid delimiter was parsed"
    end
    if r & NEWLINE > 0
        str *= ", a newline was encountered"
    end
    if r & EOF > 0
        str *= ", reached EOF"
    end
    if r & (~INVALID & INVALID_DELIMITER) > 0
        str *= ", invalid delimiter"
    end
    return str
end

codes(r) = chop(chop(string(
    ifelse(r > 0, "SUCCESS: ", "INVALID: "),
    ifelse(r & OK > 0, "OK | ", ""),
    ifelse(r & SENTINEL > 0, "SENTINEL | ", ""),
    ifelse(r & QUOTED > 0, "QUOTED | ", ""),
    ifelse(r & ESCAPED_STRING > 0, "ESCAPED_STRING | ", ""),
    ifelse(r & DELIMITED > 0, "DELIMITED | ", ""),
    ifelse(r & NEWLINE > 0, "NEWLINE | ", ""),
    ifelse(r & EOF > 0, "EOF | ", ""),
    ifelse(r & (~INVALID & INVALID_QUOTED_FIELD) > 0, "INVALID_QUOTED_FIELD | ", ""),
    ifelse(r & (~INVALID & INVALID_DELIMITER) > 0, "INVALID_DELIMITER | ", ""),
    ifelse(r & (~INVALID & OVERFLOW) > 0, "OVERFLOW | ", "")
)))

defaultdateformat(T) = Dates.dateformat""
defaultdateformat(::Type{T}) where {T <: Dates.TimeType} = Dates.default_format(T)
ptrlen(s::String) = (pointer(s), sizeof(s))


text:
function reversestring(string)
    gnirts = ""
    for i in 1:length(string)
        gnirts = gnirts * string[length(string)-i+1]
    end
    return gnirts
end
function checkpalindrome(string)
    if string == reversestring(string)
        return true
    else
        return false
    end
end
print("String to check for palindrome: ")
str = readline()
if checkpalindrome(str)
    println("$str is a palindrome!")
else
    println("$str is not a palindrome.")
end

max_stars_repo_name: ellygaytor/Projects

max_issues_repo_name: ellygaytor/Projects

max_forks_repo_name: ellygaytor/Projects


text:
# Utilities for working with HSA queues

@static if VERSION >= v"1.7-"
    """
    An HSA queue.
    Each queue is uniquely associated with an agent.
    """
    mutable struct HSAQueue
        agent::HSAAgent
        queue::Ptr{HSA.Queue}
        status::HSA.Status
        cond::Base.AsyncCondition
        @atomic active::Bool
    end
else
    """
    An HSA queue.
    Each queue is uniquely associated with an agent.
    """
    mutable struct HSAQueue
        agent::HSAAgent
        queue::Ptr{HSA.Queue}
        status::HSA.Status
        cond::Base.AsyncCondition
        active::Bool
    end
    @inline function replacefield!(value, name::Symbol, expected, desired)
        y = getfield(value, name)
        ok = y === expected
        if ok
            setfield!(value, name, desired)
        end
        return (; old=y, success=ok)
    end
end
const QUEUES = Dict{Ptr{HSA.Queue},WeakRef}()
const DEFAULT_QUEUES = IdDict{HSAAgent,HSAQueue}()

async_send(data::Ptr{Cvoid}) = ccall(:uv_async_send, Cint, (Ptr{Cvoid},), data)

function queue_error_handler(status::HSA.Status, _queue::Ptr{HSA.Queue}, queue_obj_ptr::Ptr{Cvoid})::Nothing
    if status != HSA.STATUS_SUCCESS
        queue = unsafe_pointer_to_objref(queue_obj_ptr)::HSAQueue
        queue.status = status
        async_send(queue.cond.handle)
    end

    return nothing
end

function HSAQueue(agent::HSAAgent)
    queue_size = Ref{UInt32}(0)
    getinfo(agent.agent, HSA.AGENT_INFO_QUEUE_MAX_SIZE, queue_size) |> check
    @assert queue_size[] > 0

    r_queue = Ref{Ptr{HSA.Queue}}()
    async_cond = Base.AsyncCondition()
    queue = HSAQueue(agent, Ptr{HSA.Queue}(0), HSA.STATUS_SUCCESS, async_cond, true)
    HSA.queue_create(agent.agent, queue_size[], HSA.QUEUE_TYPE_MULTI,
                     @cfunction(queue_error_handler, Cvoid, (HSA.Status, Ptr{HSA.Queue}, Ptr{Cvoid})), pointer_from_objref(queue),
                     typemax(UInt32), typemax(UInt32), r_queue) |> check
    queue.queue = r_queue[]

    lock(RT_LOCK) do
        @assert !haskey(QUEUES, queue.queue)
        QUEUES[queue.queue] = WeakRef(queue)
        _active_kernels[queue] = Vector{AMDGPU.RuntimeEvent{AMDGPU.HSAStatusSignal}}()
    end

    # Monitor queue for async errors
    # TODO: errormonitor
    queue_ptr = queue.queue
    Threads.@spawn begin
        try
            try
                wait(async_cond)
            catch err
                err isa EOFError || rethrow(err)
                return
            end
            lock(RT_LOCK) do
                if !haskey(QUEUES, queue_ptr) || QUEUES[queue_ptr].value === nothing
                    # Queue was destroyed
                    return
                end
                queue = QUEUES[queue_ptr].value
                if queue.status != HSA.STATUS_SUCCESS
                    # Queue encountered an async error, manually kill it
                    kill_queue!(queue)
                end
            end
        catch err
            @error "Queue monitor error" exception=(err,catch_backtrace())
        end
    end

    hsaref!()
    finalizer(queue) do queue
        kill_queue!(queue)
        hsaunref!()
    end
    return queue
end
HSAQueue() = HSAQueue(get_default_agent())

struct QueueError <: Exception
    queue::Ptr{HSA.Queue}
    err::HSAError
end
Base.showerror(io::IO, err::QueueError) =
    (println(io, "QueueError(Queue at $(err.queue)) due to:"); Base.showerror(io, err.err))

get_default_queue() = get_default_queue(get_default_agent())

"Gets the default queue for an agent, creating it if necessary"
function get_default_queue(agent::HSAAgent)
    lock(RT_LOCK) do
        get!(DEFAULT_QUEUES, agent) do
            HSAQueue(agent)
        end
    end
end

"Kills all kernels executing on the given queue, and destroys the queue."
function kill_queue!(queue::HSAQueue; force=false)
    if !force && !replacefield!(queue, :active, true, false).success
        return false # We didn't destroy the queue
    end
    @static if VERSION >= v"1.7-"
        @atomic queue.active = false
    else
        queue.active = false
    end

    lock(RT_LOCK) do
        delete!(QUEUES, queue.queue)
        if get(DEFAULT_QUEUES, queue.agent, nothing) == queue
            delete!(DEFAULT_QUEUES, queue.agent)
        end
        close(queue.cond)

        # Send exception to all waiter signals
        if queue.status != HSA.STATUS_SUCCESS
            err = QueueError(queue.queue, HSAError(queue.status))
            for signal in _active_kernels[queue]
                signal.exception = err
                notify(signal)
            end
        end

        delete!(_active_kernels, queue)
    end
    HSA.queue_destroy(queue.queue) |> check

    return true # We actually destroyed the queue
end


text:
# This file is a part of BAT.jl, licensed under the MIT License (MIT).

"""
    GenericDensity{F<:Function} <: AbstractDensity

**GenericDensity is deprecated and may be removed in future BAT versions**
"""
struct GenericDensity{F<:Function} <: AbstractDensity
    f::F

    @noinline function GenericDensity(f::F) where {F<:Function}
        Base.depwarn("`BAT.GenericDensity` is deprecated, use package DensityInterface.jl to turn log-density functions into BAT-compatible densities.", :GenericDensity)
        new{F}(f)
    end
end

Base.convert(::Type{GenericDensity}, f::Function) = GenericDensity(f)

@noinline function Base.convert(::Type{AbstractDensity}, f::Function)
    Base.depwarn("`convert(BAT.AbstractDensity, f::Function)` is deprecated, use `convert(AbstractDensity, logfuncdensity(g))` with a function `g` that returns the log-density value directly instead.", :convert)
    GenericDensity(f)
end

Base.parent(density::GenericDensity) = density.f

function DensityInterface.logdensityof(density::GenericDensity, v::Any)
    logvalof(density.f(v))
end

"""
    struct BAT.LFDensity{F}

*BAT-internal, not part of stable public API.*

Wraps a log-density function `log_f`.
"""
struct LFDensity{F} <: AbstractDensity
    _log_f::F
end

Base.convert(::Type{LFDensity}, density::DensityInterface.LogFuncDensity) = LFDensity(logdensityof(density))
Base.convert(::Type{AbstractDensity}, density::DensityInterface.LogFuncDensity) = convert(LFDensity, density)

@inline DensityInterface.logdensityof(density::LFDensity, x) = density._log_f(x)
@inline DensityInterface.logdensityof(density::LFDensity) = density._log_f

function Base.show(io::IO, density::LFDensity)
    print(io, Base.typename(typeof(density)).name, "(")
    show(io, density._log_f)
    print(io, ")")
end

"""
    BAT.LFDensityWithGrad{F<:Function} <: AbstractDensity

*BAT-internal, not part of stable public API.*

Constructors:

    LFDensityWithGrad(logf::Function, valgradlogf::Function)

A density defined by a function that computes it's logarithmic value at given
points, as well as a function that computes both the value and the gradient.

It must be safe to execute both functions in parallel on multiple threads and
processes.
"""
struct LFDensityWithGrad{F<:Function,G<:Function} <: AbstractDensity
    logf::F
    valgradlogf::G
end

DensityInterface.logdensityof(density::LFDensityWithGrad) = density.logf

function DensityInterface.logdensityof(density::LFDensityWithGrad, v::Any)
    density.logf(v)
end

function ChainRulesCore.rrule(::typeof(DensityInterface.logdensityof), density::LFDensityWithGrad, v)
    value, gradient = density.valgradlogf(v)
    @assert value isa Real
    function lfdwg_pullback(thunked_ΔΩ)
        ΔΩ = ChainRulesCore.unthunk(thunked_ΔΩ)
        @assert ΔΩ isa Real
        tangent = gradient * ΔΩ
        (NoTangent(), ZeroTangent(), tangent)
    end
    return value, lfdwg_pullback
end

vjp_algorithm(density::LFDensityWithGrad) = ZygoteAD()

function Base.show(io::IO, density::LFDensityWithGrad)
    print(io, Base.typename(typeof(density)).name, "(")
    show(io, density.logf)
    print(io, ")")
end