text:
function test_linalg(AT)
    @testset "linear algebra" begin
        @testset "transpose" begin
            @test compare(adjoint, AT, rand(Float32, 32, 32))
            @test compare(transpose, AT, rand(Float32, 32, 32))
        end

        @testset "permutedims" begin
            @test compare(x -> permutedims(x, (2, 1)), AT, rand(Float32, 2, 3))
            @test compare(x -> permutedims(x, (2, 1, 3)), AT, rand(Float32, 4, 5, 6))
            @test compare(x -> permutedims(x, (3, 1, 2)), AT, rand(Float32, 4, 5, 6))
        end

        @testset "issymmetric/ishermitian" begin
            n = 128
            areal = randn(n,n)/2
            aimg  = randn(n,n)/2

            @testset for eltya in (Float32, Float64, ComplexF32, ComplexF64)
                a = convert(Matrix{eltya}, eltya <: Complex ? complex.(areal, aimg) : areal)
                asym = transpose(a) + a        # symmetric indefinite
                aherm = a' + a                 # Hermitian indefinite
                @test issymmetric(asym)
                @test ishermitian(aherm)
            end
        end
    end
end


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

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

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

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

const _allowed_HardwareProfile_vmSize = ["Basic_A0", "Basic_A1", "Basic_A2", "Basic_A3", "Basic_A4", "Standard_A0", "Standard_A1", "Standard_A2", "Standard_A3", "Standard_A4", "Standard_A5", "Standard_A6", "Standard_A7", "Standard_A8", "Standard_A9", "Standard_A10", "Standard_A11", "Standard_A1_v2", "Standard_A2_v2", "Standard_A4_v2", "Standard_A8_v2", "Standard_A2m_v2", "Standard_A4m_v2", "Standard_A8m_v2", "Standard_B1s", "Standard_B1ms", "Standard_B2s", "Standard_B2ms", "Standard_B4ms", "Standard_B8ms", "Standard_D1", "Standard_D2", "Standard_D3", "Standard_D4", "Standard_D11", "Standard_D12", "Standard_D13", "Standard_D14", "Standard_D1_v2", "Standard_D2_v2", "Standard_D3_v2", "Standard_D4_v2", "Standard_D5_v2", "Standard_D2_v3", "Standard_D4_v3", "Standard_D8_v3", "Standard_D16_v3", "Standard_D32_v3", "Standard_D64_v3", "Standard_D2s_v3", "Standard_D4s_v3", "Standard_D8s_v3", "Standard_D16s_v3", "Standard_D32s_v3", "Standard_D64s_v3", "Standard_D11_v2", "Standard_D12_v2", "Standard_D13_v2", "Standard_D14_v2", "Standard_D15_v2", "Standard_DS1", "Standard_DS2", "Standard_DS3", "Standard_DS4", "Standard_DS11", "Standard_DS12", "Standard_DS13", "Standard_DS14", "Standard_DS1_v2", "Standard_DS2_v2", "Standard_DS3_v2", "Standard_DS4_v2", "Standard_DS5_v2", "Standard_DS11_v2", "Standard_DS12_v2", "Standard_DS13_v2", "Standard_DS14_v2", "Standard_DS15_v2", "Standard_DS13-4_v2", "Standard_DS13-2_v2", "Standard_DS14-8_v2", "Standard_DS14-4_v2", "Standard_E2_v3", "Standard_E4_v3", "Standard_E8_v3", "Standard_E16_v3", "Standard_E32_v3", "Standard_E64_v3", "Standard_E2s_v3", "Standard_E4s_v3", "Standard_E8s_v3", "Standard_E16s_v3", "Standard_E32s_v3", "Standard_E64s_v3", "Standard_E32-16_v3", "Standard_E32-8s_v3", "Standard_E64-32s_v3", "Standard_E64-16s_v3", "Standard_F1", "Standard_F2", "Standard_F4", "Standard_F8", "Standard_F16", "Standard_F1s", "Standard_F2s", "Standard_F4s", "Standard_F8s", "Standard_F16s", "Standard_F2s_v2", "Standard_F4s_v2", "Standard_F8s_v2", "Standard_F16s_v2", "Standard_F32s_v2", "Standard_F64s_v2", "Standard_F72s_v2", "Standard_G1", "Standard_G2", "Standard_G3", "Standard_G4", "Standard_G5", "Standard_GS1", "Standard_GS2", "Standard_GS3", "Standard_GS4", "Standard_GS5", "Standard_GS4-8", "Standard_GS4-4", "Standard_GS5-16", "Standard_GS5-8", "Standard_H8", "Standard_H16", "Standard_H8m", "Standard_H16m", "Standard_H16r", "Standard_H16mr", "Standard_L4s", "Standard_L8s", "Standard_L16s", "Standard_L32s", "Standard_M64s", "Standard_M64ms", "Standard_M128s", "Standard_M128ms", "Standard_M64-32ms", "Standard_M64-16ms", "Standard_M128-64ms", "Standard_M128-32ms", "Standard_NC6", "Standard_NC12", "Standard_NC24", "Standard_NC24r", "Standard_NC6s_v2", "Standard_NC12s_v2", "Standard_NC24s_v2", "Standard_NC24rs_v2", "Standard_NC6s_v3", "Standard_NC12s_v3", "Standard_NC24s_v3", "Standard_NC24rs_v3", "Standard_ND6s", "Standard_ND12s", "Standard_ND24s", "Standard_ND24rs", "Standard_NV6", "Standard_NV12", "Standard_NV24"]

function check_required(o::HardwareProfile)
    true
end

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

ext: jl


text:
# using LinearAlgebra
function cart2sph(x::T, y::T, z::T) where {T <: Real}
    hxy = hypot(x, y)
    r = hypot(hxy, z)
    pol = atan(hxy, z)
    az = atan(y, x)
    return [Float64(pol), Float64(az), Float64(r)]
end

function cart2sph(xyz::Point)
    x, y, z = xyz[1], xyz[2], xyz[3]
    return cart2sph(x, y, z)
end

function sph2cart(pol::T, az::T, r::T) where {T <: Real}
    x = r * sin(pol) * cos(az)
    y = r * sin(pol) * sin(az)
    z = r * cos(pol)
    return x, y, z
end
"""
Change array of vertex coordinates into array of vertex type
"""
function array_to_vertex(a::AbstractArray)
    return [Point3f0(vertex) for vertex in a]
end

"""
Change array of face into array of face type
"""
function array_to_faces(a::AbstractArray)
    return [TriangleFace{OffsetInteger{-1,UInt32}}(face) for face in a]
end

"""
Change array of array of into matrix
"""
function array_of_array_to_matrix(a::AbstractArray)
    return [a[i][j] for i = 1:length(a), j = 1:length(a[1])]
end
"""
Change matrix to array of array
"""
function matrix_to_array_of_array(mat::AbstractMatrix)
    return [mat[i, :] for i = 1:size(mat, 1)]
end # function]

"""
Change matrix to vertices
"""
function matrix_to_vertex(mat::AbstractMatrix)
    return [Point3f0(mat[i, :]) for i = 1:size(mat, 1)]
end # function

"""
Centering the coordinates of mesh vertices to their mean
"""
function centroid(v::AbstractMatrix)
    return sum(v, dims=1) / size(v, 1)
end
function centroid(verts::Vector{T1}) where {T1}
    return reduce(+, verts) / length(verts)
end

function mesh_centering(verts::Vector{T1}, faces::Vector{T2}) where {T1,T2}
    mean_coord = reduce(+, verts) / length(verts)
    new_verts = similar(verts)
    for i in eachindex(verts)
        new_verts[i] = verts[i] - mean_coord
    end
    return new_verts, faces
end
function mesh_centering(mesh::Mesh)
    verts = decompose(Point3f0, mesh)
    triangles = decompose(TriangleFace{Int}, mesh)
    verts, triangles = mesh_centering(verts, triangles)
    return Mesh(verts, faces(mesh))
end

function signed_volume_of_triangle(
    p1::T,
    p2::T,
    p3::T,
) where {T <: AbstractVector}
    return dot(p1, cross(p2, p3)) / 6.0
end

"""
These functions return the volume of mesh object
"""
function mesh_volume(
    verts::Array{Point{3,T1}},
    faces::Array{TriangleFace{T2}},
) where {T1 <: Real,T2 <: Int}
    total_vol = 0.0
    for face in faces
        vol = signed_volume_of_triangle(
            verts[face[1]],
            verts[face[2]],
            verts[face[3]],
        )
        total_vol += vol
    end
    return abs(total_vol)
end

function mesh_volume(mesh::Mesh)
    verts = decompose(Point3f0, mesh)
    faces = decompose(TriangleFace{Int}, mesh)

    return mesh_volume(verts, faces)
end
"""
These functions return the copy of mesh object with volume one
"""
function volume_normalizing(
    verts::Array{Point{3,T1}},
    faces::Array{TriangleFace{T2}},
) where {T1 <: Real,T2 <: Int}
    vol = mesh_volume(verts, faces)
    cbrt_vol = cbrt(vol)
    new_verts = similar(verts)
    for i in eachindex(verts)
        new_verts[i] = verts[i] / cbrt_vol
    end
    return new_verts, faces
end

function volume_normalizing2(
    verts::Array{Point{3,T1}},
    faces::Array{TriangleFace{T2}},
) where {T1 <: Real,T2 <: Int}
    vol = mesh_volume(verts, faces)
    cbrt_vol = cbrt(vol)
    new_verts = similar(verts)
    for i in eachindex(verts)
        new_verts[i] = verts[i] / cbrt_vol
    end
    return new_verts, faces
end

function volume_normalizing(mesh::Mesh)
    verts = decompose(Point3f0, mesh)
    triangles = decompose(TriangleFace{Int}, mesh)
    verts, triangles = volume_normalizing(verts, triangles)
    return Mesh(verts, faces(mesh))
end
"""
This function return the mesh object from verts and faces
"""
function create_mesh(verts, faces)
    faces = convert(Vector{TriangleFace{OffsetInteger{-1,UInt32}}}, faces)
    return Mesh(verts, faces)
end
"""
This function return the angle at vertex v1 from triangle(v1,v2,v3)
"""
function calculate_angle(v1, v2, v3)
    v12 = v2 - v1
    v13 = v3 - v1
    v12 = v12 / norm(v12)
    v13 = v13 / norm(v13)
    cos_val = dot(v12, v13)
    angle = acos(cos_val)
    return angle
end
# """
# This function return the internal angles of vertices
# """
# function internal_angles(verts,faces)
#     nrow,ncol = length(faces[1]),length(faces)
#     angles = Matrix{Float64}(undef,nrow,ncol)
#     for j in 1:ncol
#         f = faces[j]
#         for  i in 1:nrow
#             v1 = verts[f[i]]
#             idx2 = mod(i+1,nrow) == 0 ? 3 : mod(i+1,nrow)
#             idx3 = mod(i-1,nrow) == 0 ? 3 : mod(i-1,nrow)
#             # @show i
#             # @show idx3
#             # @show f
#             v2 = verts[f[idx2]]
#             v3 = verts[f[idx3]]
#
#             angles[i,j] = calculate_angle(v1,v2,v3)
#         end
#     end
#     return angles
# end
"""
This function return the internal angles of vertices
"""
function internal_angles_intrinsic(l)
    s23 = l[:, 1]
    s31 = l[:, 2]
    s12 = l[:, 3]

    a23 = acos.((s12.^2 + s31.^2 - s23.^2) ./ (2 * s12 .* s31))
    a31 = acos.((s23.^2 + s12.^2 - s31.^2) ./ (2 * s23 .* s12))
    a12 = acos.((s31.^2 + s23.^2 - s12.^2) ./ (2 * s31 .* s23))

    angles = [a23 a31 a12]
    return angles
end
function internal_angles(verts, faces)
    f1 = [face[1] for face in faces]
    f2 = [face[2] for face in faces]
    f3 = [face[3] for face in faces]

    s12 = norm.(verts[f2] - verts[f1])
    s13 = norm.(verts[f3] - verts[f1])
    s23 = norm.(verts[f3] - verts[f2])

    l = [s23 s13 s12]

    A = internal_angles_intrinsic(l)
end

"""
This function return the area of triangle made from vectors v1,v2, and v3)
"""
function calculate_area(v1, v2, v3)
    v12 = v2 - v1
    v13 = v3 - v1
    area = norm(cross(v12, v13)) / 2.0
    return area
end
"""
This function return the surface area)
"""
function surface_area(v::AbstractVector, f::AbstractVector)
    nv = length(v)
    nf = length(f)

    area = 0
    for i in 1:nf
        l1 = norm(v[f[i][2]] - v[f[i][3]])
        l2 = norm(v[f[i][3]] - v[f[i][1]])
        l3 = norm(v[f[i][1]] - v[f[i][2]])
        s = (l1 + l2 + l3) / 2.0
        square_area = s * (s - l1) * (s - l2) * (s - l3)
        if square_area < 0.0; square_area = 0.0 end

        area += sqrt(square_area)
    end
    return area
end

function surface_area_list(v::Matrix{<:AbstractFloat}, f::Matrix{<:Signed})
    f1 = @view f[:, 1]
    f2 = @view f[:, 2]
    f3 = @view f[:, 3]

    l1 = sqrt.(sum((view(v, f2, :) - view(v, f3, :)).^2, dims=2))
    l2 = sqrt.(sum((view(v, f3, :) - view(v, f1, :)).^2, dims=2))
    l3 = sqrt.(sum((view(v, f1, :) - view(v, f2, :)).^2, dims=2))
    s = (l1 + l2 + l3) / 2.0
    squared_area = s .* (s - l1) .* (s - l2) .* (s - l3)
    area = sqrt.(squared_area)
    return area
end
function surface_area(v::Matrix{<:AbstractFloat}, f::Matrix{<:Signed})
    area = surface_area_list(v, f)
    return sum(area)
end
function surface_area2(v::AbstractVector, f::AbstractVector)
    f1 = [face[1] for face in f]
    f2 = [face[2] for face in f]
    f3 = [face[3] for face in f]
    # F x 3 matrices of triangle edge vectors, named after opposite vertices
    v1 = [v[f3[i]] - v[f2[i]] for i in eachindex(f)]
    v2 = [v[f1[i]] - v[f3[i]] for i in eachindex(f)]
    # v3 = [v[f2[i]] - v[f1[i]] for i in eachindex(f)];
    n = cross.(v1, v2)
    area = norm.(n) / 2.0
    return sum(area)
end
function surface_area2(v::Matrix{<:AbstractFloat}, f::Matrix{<:Signed})
    f1 = @view f[:, 1]
    f2 = @view f[:, 2]
    f3 = @view f[:, 3]
    # F x 3 matrices of triangle edge vectors, named after opposite vertices
    v1 = view(v, f3, :) - view(v, f2, :)
    v2 = view(v, f1, :) - view(v, f3, :)
    # v3 = v[f2,:] - v[f1,:];
    n = [cross(v1[i, :], v2[i, :]) for i = 1:size(v1, 1)]
    area = norm.(n) / 2.0
    return sum(area)
end

function mesh_reorientation(mesh::Mesh, T, N, B)
    verts = decompose(Point3f0, mesh)
    mat_verts = array_of_array_to_matrix(verts)
    new_mat_verts = [T N B] \ mat_verts'
    new_mat_verts = new_mat_verts'
    new_verts = matrix_to_vertex(new_mat_verts)
    return Mesh(new_verts, faces(mesh))
end
function mesh_reorientation(mat_verts::AbstractMatrix, T, N, B)
    new_mat_verts = [T N B] \ mat_verts'
    return new_mat_verts'
end


text:
using GPs
using Base.Test
using StaticArrays

#--- Metric tests

@testset "Utils tests" begin

    for x in [toSVector([1, 2, 3]),
              toSVector([SVector(1), SVector(2), SVector(3)]),
              toSVector([1 2 3]),
              toSVector(1:3)]
        @test x == [SVector(1), SVector(2), SVector(3)]
        @test isa(x, Vector{SVector{1, Int}})
    end
    for x in [toMatrix([1, 2, 3]),
              toMatrix([SVector(1), SVector(2), SVector(3)]),
              toMatrix([1 2 3])]
        @test x == [1 2 3]
        @test isa(x, Matrix)
    end

end

@testset "Metric tests" begin

    @testset "Isotropic Euclidean" begin

        x1 = @SVector [1.3, 2.2, 4.2]
        x2 = @SVector [-0.2, 2.3, 2.1]
        iso = GPs.Metrics.Euclidean(0.1)

        @test GPs.hyperparametercount(iso) == 1
        GPs.sethyperparameters!(iso, 0.5)
        @test GPs.gethyperparameters(iso) == [0.5]

        @test GPs.Metrics.sqrdistance(iso, x1, x2) ≈ 26.68
        ans1 = GPs.Metrics.sqrdistance_dθ(iso, x1, x2)
        @test ans1 ≈ [-106.72] # check result
        @test isa(ans1, SVector) # check that types have propogated correctly
        ans2 = GPs.Metrics.sqrdistance_dx1(iso, x1, x2)
        @test ans2 ≈ [12.00, -0.80, 16.80]
        @test isa(ans2, SVector)
        ans3 = GPs.Metrics.sqrdistance_dx2(iso, x1, x2)
        @test ans3 ≈ [-12.00, 0.80, -16.80]
        @test isa(ans3, SVector)

    end

    @testset "Anisotropic Euclidean" begin

        x1 = @SVector [1.3, 2.2, 4.2]
        x2 = @SVector [-0.2, 2.3, 2.1]
        aniso = GPs.Metrics.Euclidean(MVector(0.5, 0.2, 1.2))

        @test GPs.hyperparametercount(aniso) == 3
        GPs.sethyperparameters!(aniso, [0.5, 0.2, 1.3])
        @test GPs.gethyperparameters(aniso) == MVector(0.5, 0.2, 1.3)

        @test GPs.Metrics.sqrdistance(aniso, x1, x2) ≈ 11.85946746
        ans1 = GPs.Metrics.sqrdistance_dθ(aniso, x1, x2)
        @test ans1 ≈ [-36.0, -2.5, -4.014565316]
        @test isa(ans1, SVector)
        ans2 = GPs.Metrics.sqrdistance_dx1(aniso, x1, x2)
        @test ans2 ≈ [12.00, -5.00, 2.485207100]
        @test isa(ans2, SVector)
        ans3 = GPs.Metrics.sqrdistance_dx2(aniso, x1, x2)
        @test ans3 ≈ [-12.00, 5.00, -2.485207100]
        @test isa(ans3, SVector)

    end

end

@testset "Kernel tests" begin

    @testset "Squared-exponential kernel" begin

        iso = GPs.Metrics.Euclidean(0.5)
        sqrexp = GPs.Kernels.SqrExponential(iso, 0.125)

        x1 = @SVector [1.3, 2.2, 4.2]
        x2 = @SVector [-0.2, 2.3, 2.1]

        covar = GPs.Kernels.covariance(sqrexp, x1, x2)
        @test covar ≈ 2.5138053e-08

        ans1 = GPs.Kernels.covariance_dθ(sqrexp, x1, x2)
        @test ans1 ≈ [4.022088479e-07, 1.341366507e-06]
        @test isa(ans1, SVector)
        ans2 = GPs.Kernels.covariance_dx1(sqrexp, x1, x2)
        @test ans2 ≈ [-1.508283182e-07, 1.005522113e-08, -2.111596459e-07]
        @test isa(ans2, SVector)
        ans3 = GPs.Kernels.covariance_dx2(sqrexp, x1, x2)
        @test ans3 ≈ [1.508283182e-07, -1.005522113e-08, 2.111596459e-07]

        @test GPs.hyperparametercount(sqrexp) == 2
        GPs.sethyperparameters!(sqrexp, [0.55, 1.3])
        @test GPs.gethyperparameters(sqrexp) == [0.55, 1.3]

        @test GPs.Kernels.covariance(sqrexp, x1, x2) ≈ 0.042043656
    end

end

@testset "Covariance matrix tests" begin

    iso = GPs.Metrics.Euclidean(0.5)
    sqrexp = GPs.Kernels.SqrExponential(iso, 0.125)

    x1 = [SVector(1.3, 2.2, 4.2), SVector(-0.2, 1.2, -0.1)]
    x2 = [SVector(-0.2, 2.3, 2.1)]

    @test GPs.Kernels.covariance_matrix(sqrexp, x1, x2) ≈ [2.513805299e-08; 8.68673944e-08]
    M1 = GPs.Kernels.covariance_matrix(sqrexp, x1)
    @test M1 ≈ [0.015625 2.045012434e-21; 2.045012434e-21 0.015625]

    M2 = zeros(2, 2)
    GPs.Kernels.covariance_matrix!(UpperTriangular(M2), sqrexp, x1)
    M1[2, 1] = 0.0
    @test M1 == M2

end

@testset "Gaussian Process (Gaussian likelihood) tests" begin

    n = 21
    x = linspace(0, 1, n)
    y = sin.(2π*x)

    σn = 0.3
    σf = 0.1
    ℓ = 0.2

    iso = GPs.Metrics.Euclidean(ℓ)
    sqrexp = GPs.Kernels.SqrExponential(iso, σf)
    training = GPs.DataSets.TrainingSet(x, y)
    gp = GPs.GaussianLikelihoods.GP(sqrexp, σn, training)

    @test isa(gethyperparameters(gp), SVector{3, Float64})

    # Calculate the mean between the training points
    xx = 0.5*(x[2:end] + x[1:end-1])

    # From Matlab GPML toolbox (ell = log(0.2), sf = log(0.1), lik_sn = log(0.3))
    #   n = 21;
    #   x = linspace(0, 1, n);
    #   y = sin(x);
    #   xx = 0.5*(x(2:end) + x(1:end-1))
    #   meanfunc = [];                    % empty: don't use a mean function
    #   covfunc = @covSEiso;              % Squared Exponental covariance function
    #   likfunc = @likGauss;              % Gaussian likelihood
    #   hyp = struct('mean', [], 'cov', [log(0.2) log(0.1)], 'lik', log(0.1));
    #   [mu s2] = gp(hyp, @infGaussLik, meanfunc, covfunc, likfunc, x, y, xx);

    μ₀ = [0.202764557432603 0.2532874969517 0.298316602533046 0.330784847496406 0.34404680063947 0.333129086919798 0.295761669144916 0.232966591637241 0.149093996664424 0.0513161114615125 -0.0513161114615124 -0.149093996664424 -0.232966591637241 -0.295761669144916 -0.333129086919798 -0.34404680063947 -0.330784847496406 -0.298316602533047 -0.2532874969517 -0.202764557432603]
    σ₀² = [0.0968914136106435,0.0964262935050147,0.0961104866357293,0.0959299264929463,0.0958451743017393,0.0958140514271935,0.0958060359999692,0.0958050354333621,0.0958051126669153,0.0958050947257113,0.0958050947257113,0.0958051126669153,0.0958050354333621,0.0958060359999692,0.0958140514271935,0.0958451743017393,0.0959299264929463,0.0961104866357293,0.0964262935050147,0.0968914136106435]

    μ = mean(gp, xx)
    @test μ ≈ μ₀
    σ² = var(gp, xx)
    @test σ² ≈ σ₀²
    K = cov(gp, xx)
    @test diag(K) ≈ σ₀²

    # Test multiple outputs
    training2 = GPs.DataSets.TrainingSet(x, [y y]')
    gp2 = GPs.GaussianLikelihoods.GP(sqrexp, σn, training2)

    μ2 = mean(gp2, xx)
    @test μ2[1, :] == μ2[2, :]
    @test μ2[1:1, :] ≈ μ₀

    # Test likelihood calculations
    x = [0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95, 1.0]
    y = [0.107533 0.675794 0.136016 0.981452 1.01481 0.738462 0.864338 0.877542 1.30346 0.862904 -0.269977 0.297968 -0.442704 -0.821628 -0.808108 -1.04099 -0.975885 -0.511077 -0.305978 -0.0255785 0.134299; -0.241497 0.452465 0.913832 0.906796 1.158 1.14538 0.890368 0.867791 0.430329 0.486696 -0.229414 -0.522791 -0.749685 -1.39787 -0.66338 -0.934962 -1.10204 -0.534957 -0.930089 -0.329465 -0.0482894]
    σn = 0.1
    σf = 0.1
    ℓ = 0.2

    iso = GPs.Metrics.Euclidean(ℓ)
    sqrexp = GPs.Kernels.SqrExponential(iso, σf)
    GP = GPs.GaussianLikelihoods.GP
    trainingy = GPs.DataSets.TrainingSet(x, y)
    trainingy1 = GPs.DataSets.TrainingSet(x, y[1, :])
    trainingy2 = GPs.DataSets.TrainingSet(x, y[2, :])

    @test round(loglikelihood(GP(sqrexp, σn, trainingy1)), 2) ≈ round(136.272815830685, 2)
    @test round(loglikelihood(GP(sqrexp, σn, trainingy2)), 2) ≈ round(132.426903505274, 2)
    @test round(loglikelihood(GP(sqrexp, σn, trainingy)), 2) ≈ round(mean([136.272815830685, 132.426903505274]), 2)

    @test loglikelihood_dθ(GP(sqrexp, σn, trainingy)) ≈ [-1452.713883, -1531.4908, 306.8887837]

    # Test number type propogation
    n = 21
    x = convert(Vector{Float32}, linspace(0, 1, n))
    y = sin.(Float32(2)*π*x)

    σn = Float32(0.3)
    σf = Float32(0.1)
    ℓ = Float32(0.2)

    iso = GPs.Metrics.Euclidean(ℓ)
    sqrexp = GPs.Kernels.SqrExponential(iso, σf)
    training = GPs.DataSets.TrainingSet(x, y)
    gp = GPs.GaussianLikelihoods.GP(sqrexp, σn, training)

    # Calculate the mean between the training points
    xx = Float32(0.5)*(x[2:end] + x[1:end-1])

    μ = mean(gp, xx)
    @test isa(μ, Matrix{Float32})
    σ² = var(gp, xx)
    @test isa(σ², Vector{Float32})
    K = cov(gp, xx)
    @test isa(K, Matrix{Float32})
    L = loglikelihood(gp)
    @test isa(L, Float32)

end


text:
module TestingLinearInterpolatorsInterpolations

using TwoDimensional

using LazyAlgebra
using LazyAlgebra: Adjoint, Direct

using LinearInterpolators
using LinearInterpolators: getcoefs, limits

using Printf
using Test
using Statistics: mean

@testset "Unidimensional" begin
    ker = @inferred LinearSpline(Float32)
    src = [0, 0, 1, 0, 0] # Vector{Int}
    x = LinRange(-1,7,81) # Float64
    out = @inferred apply(ker, x, src) # different eltypes
    @test out[x .== 2.5][1] == 0.5

    out2 = similar(out)
    apply!(out2, ker, x, src)
    @test out2 == out

    apply!(2, Direct, ker, x, src, 0, out2) # α=2 β=0
    @test out2 == 2*out

    out2[:] .= 1
    @inferred apply!(1, Direct, ker, x, src, 2, out2) # α=1 β=2
    @test out2 == out .+ 2

    out2[:] .= 3
    @inferred apply!(0, Direct, ker, x, src, 2, out2) # α=0 β=2
    @test all(==(6), out2)

    adj = @inferred apply(Adjoint, ker, x, out, length(src))
    adj2 = similar(adj)
    @inferred apply!(2, Adjoint, ker, x, out, 0, adj2) # α=2 β=0
    @test adj2 == 2*adj

    adj2[:] .= 1
    @inferred apply!(1, Adjoint, ker, x, out, 2, adj2) # α=1 β=2
    @test adj2 ≈ adj .+ 2

    adj2[:] .= 3
    @inferred apply!(0, Adjoint, ker, x, out, 2, adj2) # α=0 β=2
    @test all(==(6), adj2)
end

distance(a::Real, b::Real) = abs(a - b)

distance(a::NTuple{2,Real}, b::NTuple{2,Real}) =
    hypot(a[1] - b[1], a[2] - b[2])

distance(A::AffineTransform2D, B::AffineTransform2D) =
    max(abs(A.xx - B.xx), abs(A.xy - B.xy), abs(A.x - B.x),
        abs(A.yx - B.yx), abs(A.yy - B.yy), abs(A.y - B.y))

distance(A::AbstractArray{Ta,N}, B::AbstractArray{Tb,N}) where {Ta,Tb,N} =
    mean(abs.(A - B))

shortname(::Nothing) = ""
shortname(m::RegexMatch) = m.captures[1]
shortname(::Type{T}) where {T} = shortname(string(T))
shortname(str::AbstractString) =
    shortname(match(r"([_A-Za-z][_A-Za-z0-9]*)([({]|$)", str))

struct CustomCoordinate{T<:Real}
    val::T
end

LinearInterpolators.convert_coordinate(T::Type, c::CustomCoordinate) =
    convert(T, c.val)

kernels = (RectangularSpline(), LinearSpline(), QuadraticSpline(),
           CubicSpline(), CatmullRomSpline(), KeysSpline(-0.4),
           MitchellNetravaliSpline(), LanczosKernel(4), LanczosKernel(6))

conditions = (Flat, SafeFlat)

const VERBOSE = true
sub = 3
x = range(-1, stop=1, length=100*sub + 1);
xsub = x[1:sub:end];
y = cos.(2 .* x .* (x .+ 2));
ysub = y[1:sub:end];
t = range(1, stop=length(xsub), length=length(x));

@testset "`getcoefs` method" begin
    T = Float64
    ker = CatmullRomSpline(T)
    len = 4
    lim = limits(ker, len)
    for x in (0, .1, -1.2)
        @test getcoefs(ker, lim, x) == getcoefs(ker, lim, CustomCoordinate(x))
    end
end

@testset "TabulatedInterpolators" begin
    tol = 1e-14

    for K in kernels,
        C in conditions
        tol = isa(K, RectangularSpline) ? 0.02 : 0.006
        ker = C(K)
        T = @inferred TabulatedInterpolator(ker, t, length(xsub))
        @inferred T(ysub)
        @test distance(T(ysub), T*ysub) ≤ 0
        err = distance(T*ysub, y)
        if VERBOSE
            print(shortname(typeof(K)),"/",shortname(C)," max. err = ")
            @printf("%.3g\n", err)
        end
        @test err ≤ tol
        @test distance(T*ysub, apply(ker,t,ysub)) ≤ 1e-15

    end
end

@testset "SparseInterpolators" begin
    for K in kernels, C in conditions
        tol = isa(K, RectangularSpline) ? 0.02 : 0.006
        ker = C(K)
        S = @inferred SparseInterpolator(ker, t, length(xsub))
        T = @inferred TabulatedInterpolator(ker, t, length(xsub))
        @inferred S(ysub)
        @inferred T(ysub)
        @test distance(S(ysub), S*ysub) ≤ 1e-15
        @test distance(S(ysub), T(ysub)) ≤ 1e-15
        @test distance(S(ysub), y) ≤ tol
        @test distance(S(ysub), apply(ker,t,ysub)) ≤ 1e-15
    end
    let ker = kernels[1], T = Float32
        @test_deprecated SparseInterpolator(T, ker, t, length(xsub))
    end
end

@testset "TwoDimensionalTransformInterpolator" begin
    kerlist = (RectangularSpline(), LinearSpline(), CatmullRomSpline())
    rows = (11,16)
    cols = (10,15)

    # Define a rotation around c.
    c = (5.6, 6.4)
    R = c + rotate(AffineTransform2D{Float64}() - c, 0.2)

    for ker1 in kerlist, ker2 in kerlist
        A = @inferred TwoDimensionalTransformInterpolator(rows, cols, ker1, ker2, R)
        x = randn(cols)
        y = randn(rows)
        #x0 = randn(cols)
        #y0 = Array{Float64}(undef, rows)
        #x1 = Array{Float64}(undef, cols)
        @test vdot(A*x, y) ≈ vdot(x, A'*y)
    end

end

end # module


text:
# Matrix Product States
# #+HTML: Source
# #+HTML:

# [[file:~/.julia/dev/MatrixProductStates/README.org::*Matrix%20Product%20States][Matrix Product States:1]]
"""
    MPS{L, T<:Number}

Matrix product state on L sites.

The `i`th tensor in the state has indices `[aⁱ⁻¹, aⁱ, σⁱ]` where
`(aⁱ⁻¹, aⁱ)` are bond indices and `σⁱ` is the physical index.

A four site MPS would be diagrammatically represented

    σ¹          σ²          σ³          σ⁴
    |           |           |           |
    •--(a¹ a¹)--•--(a² a²)--•--(a³ a³)--•

Note that `a⁰` and `aᴸ` must be of dimension 1.
"""
struct MPS{L, T<:Number}
    tensors::Vector{Array{T,3}}
end

Base.isequal(ψ::MPS, ϕ::MPS)     = (isequal(ψ.tensors, ϕ.tensors))
Base.isapprox(ψ::MPS, ϕ::MPS)   = isapprox(ψ.tensors, ϕ.tensors)

Base.eltype(::Type{MPS{L, T}}) where {L, T} = T

Base.length(::MPS{L, T}) where {L, T} = L

Base.size(::MPS{L, T}) where {L, T} = (L,)
Base.getindex(ψ::MPS, i::Int) = getindex(ψ.tensors, i)

Base.:(*)(ψ::MPS{L, T}, x::Number) where {L, T} = MPS{L,T}(ψ.tensors .* x)
Base.:(*)(x::Number, ψ::MPS) = ψ * x
Base.:(/)(ψ::MPS{L,T}, x::Number) where {L, T} = MPS{L,T}(ψ.tensors ./ x)
Base.copy(ψ::MPS{L, T}) where {L, T} = MPS{L,T}(copy.(ψ.tensors))

function Base.randn(::Type{MPS{L, T}}, D::Int, d::Int) where {L, T}
    tensors = [randn(1, D, d), [randn(D, D, d) for _ in 2:(L-1)]..., randn(D, 1, d)]
    MPS{L, T}(tensors) |> leftcanonical |> rightcanonical
end

"""
    MPS(vs::Vector{Vector})
Create an `MPS` representing a product state (all bonds have dimension 1),
where each site is described by the corresponding element of `vs`.
"""
function MPS(vs::Vector{Vector{T}}) where {T}
    L = length(vs)

    tensrs = Vector{Array{T,3}}(undef, L)
    for i in 1:L
        tensrs[i] = reshape(copy(vs[i]), 1, 1, :)
    end

    MPS{L,T}(tensrs)
end

"""
    MPS(v::Vector, L)
Create an `MPS` for `L` sites representing a uniform product state (all bonds
have dimension 1), where each site is described by `v`.
"""
MPS(v::Vector, L) = MPS([v for _ in 1:L])

function Base.show(io::IO, ::MIME"text/plain", ψ::MPS{L, T}) where {L, T}
    d = length(ψ.tensors[2][1, 1, :])
    bonddims = [size(ψ[i][:, :, 1]) for i in 1:L]
    println(io, "Matrix product state on $L sites")
    _show_mps_dims(io, L, d, bonddims)
end

function Base.show(ψ::MPS{L, T}) where {L, T}
    d = length(ψ.tensors[2][1, 1, :])
    bonddims = [size(ψ[i][:, :, 1]) for i in 1:L]
    println("Matrix product state on $L sites")
    _show_mps_dims(L, d, bonddims)
end

function _show_mps_dims(io::IO, L, d, bonddims)
    println(io, "  Physical dimension: $d")
    print(io, "  Bond dimensions:   ")
    if L > 8
        for i in 1:8
            print(io, bonddims[i], " × ")
        end
        print(io, " ... × ", bonddims[L])
    else
        for i in 1:(L-1)
            print(io, bonddims[i], " × ")
        end
        print(io, bonddims[L])
    end
end

function Base.show(io::IO, ψ::MPS{L, T}) where {L, T}
    print(io, "MPS on $L sites")
end
# Matrix Product States:1 ends here

# #+HTML: Adjoint MPS
# #+HTML:

# [[file:~/.julia/dev/MatrixProductStates/README.org::*Matrix%20Product%20States][Matrix Product States:2]]
function Base.adjoint(ψ::MPS{L, T}) where {L,T}
    Adjoint{T, MPS{L, T}}(ψ)
end

function Base.show(io::IO, ::MIME"text/plain", ψ::Adjoint{T, MPS{L, T}}) where {L, T}
    d = length(ψ.parent[2][1, 1, :])
    bonddims = reverse([reverse(size(ψ.parent[i][:, :, 1])) for i in 1:L])
    println(io, "Adjoint matrix product state on $L sites")
    _show_mps_dims(io, L, d, bonddims)
end

function Base.show(io::IO, ψ::Adjoint{T, MPS{L, T}}) where {L, T}
    print(io, "Adjoint MPO on $L sites")t
end

Base.size(::Adjoint{T, MPS{L, T}}) where {L, T} = (1, L)

function Base.getindex(ψ::Adjoint{T, MPS{L, T}}, args...) where {L, T}
    out = getindex(reverse(ψ.parent.tensors), args...)
    permutedims(conj.(out), (2, 1, 3))
end

adjoint_tensors(ψ::MPS) = reverse(conj.(permutedims.(ψ.tensors, [(2, 1, 3)])))
# Matrix Product States:2 ends here


text:
const TEST_NAME = "MT19937"

include("common.jl")

using RandomNumbers.MersenneTwisters

r = MT19937(123)

test_all(r, 100_000_000)


text:
using JuMP, Random

Random.seed!(1234);

a = rand(1:100, 5, 5)
c = rand(1:100, 5, 5)
b = rand(1:100, 5)
d = rand(1:100, 5)

model = Model()
@variable(model, x[1:5])
@variable(model, y, Bin)
@constraint(model, a * x .>= y .* b)
@constraint(model, c * x .>= (1 - y) .* d);

a = rand(1:100, 5, 5)
b = rand(1:100, 5)
m = rand(10000:11000, 5)

model = Model()
@variable(model, x[1:5])
@variable(model, z, Bin)
@constraint(model, a * x .<=  b .+ (m .* (1 - z)));
# If z was a regular Julia variable, we would not have had to use the vectorized dot operator

model = Model()

@variable(model, x)
@variable(model, y)
@constraint(model, x in MOI.ZeroOne())
@constraint(model, y in MOI.Integer());

l = 7.45
u = 22.22
@variable(model, a)
@constraint(model, a in MOI.Semicontinuous(l, u))

l = 5
u = 34
@variable(model, b)
@constraint(model, b in MOI.Semiinteger(l, u))

@variable(model, u[1:3])
@constraint(model, u in MOI.SOS1([1.0, 2.0, 3.0]))

@variable(model, v[1:3])
@constraint(model, v in MOI.SOS2([3.0, 1.0, 2.0]))


text:
using HttpCommon
using Base.Test
using Compat
import Compat: UTF8String
# headers
@test isa(HttpCommon.headers(), Headers)

# Request
@test sprint(show, Request()) == "Request(:/, 0 headers, 0 bytes in body)"
@test sprint(show, Request("GET", "/get", HttpCommon.headers(), "")) ==
        "Request(:/, 4 headers, 0 bytes in body)"

# Cookie
@test sprint(show, Cookie("test","it")) ==
        "Cookie(test, it, 0 attributes)"

# Response
h = HttpCommon.headers()
@test sprint(show, Response(200, h, "HttpCommon")) ==
        "Response(200 OK, 4 headers, 10 bytes in body)"
@test sprint(show, Response(200, h, UInt8[1,2,3])) ==
        "Response(200 OK, 4 headers, 3 bytes in body)"
@test sprint(show, Response(200, h)) ==
        "Response(200 OK, 4 headers, 0 bytes in body)"
@test sprint(show, Response(200, "HttpCommon")) ==
        "Response(200 OK, 4 headers, 10 bytes in body)"
@test sprint(show, Response(200, UInt8[1,2,3])) ==
        "Response(200 OK, 4 headers, 3 bytes in body)"
@test sprint(show, Response("HttpCommon", h)) ==
        "Response(200 OK, 4 headers, 10 bytes in body)"
@test sprint(show, Response(UInt8[1,2,3], h)) ==
        "Response(200 OK, 4 headers, 3 bytes in body)"
@test sprint(show, Response("HttpCommon")) ==
        "Response(200 OK, 4 headers, 10 bytes in body)"
@test sprint(show, Response(UInt8[1,2,3])) ==
        "Response(200 OK, 4 headers, 3 bytes in body)"
@test sprint(show, Response(200)) ==
        "Response(200 OK, 4 headers, 0 bytes in body)"
@test sprint(show, Response()) ==
        "Response(200 OK, 4 headers, 0 bytes in body)"

# Escape HTML
@test escapeHTML(" foo bar") ==
        "<script type='text/javascript'>alert('sucker');</script> foo bar"

# Parse URL query strings
@test parsequerystring("foo=%3Ca%20href%3D%27foo%27%3Ebar%3C%2Fa%3Erun%26%2B%2B&bar=123") ==
        Dict("foo" => "barrun&++", "bar" => "123")

begin
  substrings = split(UTF8String("a%20=1&b=%202,b,c"), ",")
  @test parsequerystring(substrings[1]) == Dict("a " => "1", "b" => " 2")
end
@test parsequerystring("") == Dict()
@test_throws ArgumentError parsequerystring("looknopairs")
@test_throws ArgumentError parsequerystring("good=pair&badpair")


text:
module ChaoticCommunications

using DifferentialEquations
using Statistics
using DocStringExtensions

include("maps.jl")
include("blocks.jl")

end # module


text:
module REPL

    # compatibility for julia version 0.2
    if !isdefined(:put!)
        const put! = put
        const take! = take
    end

    export StreamREPL, BasicREPL

    abstract AbstractREPL

    type REPLBackend
        repl_channel::RemoteRef
        response_channel::RemoteRef
        ans
    end

    using Base.Meta

    function eval_user_input(ast::ANY, backend)
        iserr, lasterr, bt = false, (), nothing
        while true
            try
                if iserr
                    put!(backend.response_channel,(lasterr,bt))
                    iserr, lasterr = false, ()
                else
                    ast = expand(ast)
                    ans = Base.Meta.quot(backend.ans)
                    eval(Main,:(ans=$(ans)))
                    value = eval(Main,ast)
                    backend.ans = value
                    put!(backend.response_channel,(value,nothing))
                end
                break
            catch err
                if iserr
                    println("SYSTEM ERROR: Failed to report error to REPL frontend")
                    println(err)
                end
                iserr, lasterr = true, err
                bt = catch_backtrace()
            end
        end
    end

    function parse_input_line(s::String)
        # s = bytestring(s)
        # (expr, pos) = parse(s, 1)
        # (ex, pos) = ccall(:jl_parse_string, Any,
        #                   (Ptr{Uint8},Int32,Int32),
        #                   s, int32(pos)-1, 1)
        # if !is(ex,())
        #     throw(ParseError("extra input after end of expression"))
        # end
        # expr
        ccall(:jl_parse_input_line, Any, (Ptr{Uint8},), s)
    end

    function start_repl_backend(repl_channel, response_channel)
        backend = REPLBackend(repl_channel,response_channel,nothing)
        @async begin
            # include looks at this to determine the relative include path
            # nothing means cwd
            while true
                tls = task_local_storage()
                tls[:SOURCE_PATH] = nothing
                (ast,show_value) = take!(backend.repl_channel)
                if show_value == -1
                    # exit flag
                    break
                end
                eval_user_input(ast, backend)
            end

        end
    end

    function display_error(io::IO, er, bt)
        Base.with_output_color(:red, io) do io
            print(io, "ERROR: ")
            Base.showerror(io, er, bt)
        end
    end

    import Base: Display, display, writemime

    immutable REPLDisplay <: Display
        repl::AbstractREPL
    end
    function display(d::REPLDisplay, ::MIME"text/plain", x)
        io = outstream(d.repl)
        write(io,answer_color(d.repl))
        writemime(io, MIME("text/plain"), x)
        println(io)
    end
    display(d::REPLDisplay, x) = display(d, MIME("text/plain"), x)

    function print_response(d::REPLDisplay,errio::IO,r::AbstractREPL,val::ANY, bt, show_value, have_color)
        while true
            try
                if !is(bt,nothing)
                    display_error(errio,val,bt)
                    println(errio)
                    iserr, lasterr = false, ()
                else
                    if !is(val,nothing) && show_value
                        try display(d,val)
                        catch err
                            println(errio,"Error showing value of type ", typeof(val), ":")
                            rethrow(err)
                        end
                    end
                end
                break
            catch err
                if !is(bt,nothing)
                    println(errio,"SYSTEM: show(lasterr) caused an error")
                    break
                end
                val = err
                bt = catch_backtrace()
            end
        end
    end

    using Terminals
    using LineEdit

    import LineEdit: char_move_left, char_move_word_left, CompletionProvider, completeLine

    type LineEditREPL <: AbstractREPL
        t::TextTerminal
        prompt_color::String
        input_color::String
        answer_color::String
        shell_color::String
        help_color::String
        in_shell::Bool
        in_help::Bool
        consecutive_returns
    end
    outstream(r::LineEditREPL) = r.t

    LineEditREPL(t::TextTerminal) =  LineEditREPL(t,julia_green,Base.input_color(),Base.answer_color(),Base.text_colors[:red],Base.text_colors[:yellow],false,false,0)

    type REPLCompletionProvider <: CompletionProvider
        r::LineEditREPL
    end

    type ShellCompletionProvider <: CompletionProvider
        r::LineEditREPL
    end

    using REPLCompletions

    function completeLine(c::REPLCompletionProvider,s)
        partial = bytestring(s.input_buffer.data[1:position(s.input_buffer)])
        ret, range = completions(partial,endof(partial))
        return (ret,partial[range])
    end

    function completeLine(c::ShellCompletionProvider,s)
        # First parse everything up to the current position
        partial = bytestring(s.input_buffer.data[1:position(s.input_buffer)])
        ret, range = shell_completions(partial,endof(partial))
        return (ret, partial[range])
    end

    import LineEdit: HistoryProvider, add_history, history_prev, history_next, history_search

    type REPLHistoryProvider <: HistoryProvider
        history::Array{String,1}
        history_file
        cur_idx::Int
        last_buffer::IOBuffer
        last_mode
        mode_mapping
        modes::Array{Uint8,1}
    end

    function hist_from_file(file,mode_mapping)
        hp = REPLHistoryProvider(String[],file,0,IOBuffer(),nothing,mode_mapping,Uint8[])
        seek(file,0)
        while !eof(file)
            b = readuntil(file,'\0')
            if uint8(b[1]) in keys(mode_mapping)
                push!(hp.modes,uint8(b[1]))
                push!(hp.history,b[2:(end-1)]) # Strip trailing \0
            else # For history backward compatibility
                push!(hp.modes,0)
                push!(hp.history,b[1:(end-1)]) # Strip trailing \0
            end
        end
        seekend(file)
        hp
    end

    function mode_idx(hist::REPLHistoryProvider,mode)
        c::Uint8 = 0
        for (k,v) in hist.mode_mapping
            if k == uint8('\0')
                continue
            elseif v == mode
                c = k
                break
            end
        end
        @assert c != 0
        return c
    end

    function add_history(hist::REPLHistoryProvider,s)
        # bytestring copies
        str = bytestring(pointer(s.input_buffer.data),s.input_buffer.size)
        if isempty(strip(str)) # Do not add empty strings to the history
            return
        end
        push!(hist.history,str)
        c = mode_idx(hist,LineEdit.mode(s))
        push!(hist.modes,c)
        write(hist.history_file,c)
        write(hist.history_file,str)
        write(hist.history_file,'\0')
        flush(hist.history_file)
    end

    function history_adjust(hist::REPLHistoryProvider,s)
        if 0 < hist.cur_idx <= length(hist.history)
            hist.history[hist.cur_idx] = LineEdit.input_string(s)
            hist.modes[hist.cur_idx] = mode_idx(hist,LineEdit.mode(s))
        end
    end

    function history_prev(s::LineEdit.MIState,hist::REPLHistoryProvider)
        if hist.cur_idx > 1
            if hist.cur_idx == length(hist.history)+1
                hist.last_mode = LineEdit.mode(s)
                hist.last_buffer = copy(LineEdit.buffer(s))
            else
                history_adjust(hist,s)
            end
            hist.cur_idx-=1
            LineEdit.transition(s,hist.mode_mapping[hist.modes[hist.cur_idx]])
            LineEdit.replace_line(s,hist.history[hist.cur_idx])
            LineEdit.refresh_line(s)
        else
            Terminals.beep(LineEdit.terminal(s))
        end
    end

    function history_next(s::LineEdit.MIState,hist::REPLHistoryProvider)
        if hist.cur_idx < length(hist.history)
            history_adjust(hist,s)
            hist.cur_idx+=1
            LineEdit.transition(s,hist.mode_mapping[hist.modes[hist.cur_idx]])
            LineEdit.replace_line(s,hist.history[hist.cur_idx])
            LineEdit.refresh_line(s)
        elseif hist.cur_idx == length(hist.history)
            hist.cur_idx+=1
            buf = hist.last_buffer
            hist.last_buffer = IOBuffer()
            LineEdit.transition(s,hist.last_mode)
            LineEdit.replace_line(s,buf)
            LineEdit.refresh_line(s)
        else
            Terminals.beep(LineEdit.terminal(s))
        end
    end

    function history_search(hist::REPLHistoryProvider,query_buffer::IOBuffer,response_buffer::IOBuffer,backwards::Bool=false, skip_current::Bool=false)
        if !(query_buffer.ptr > 1)
            #truncate(response_buffer,0)
            return true
        end

        # Alright, first try to see if the current match still works
        searchdata = bytestring(query_buffer.data[1:(query_buffer.ptr-1)])
        pos = position(response_buffer)
        if !skip_current && !((response_buffer.size == 0) || (pos+query_buffer.ptr-2 == 0)) &&
            (response_buffer.size >= (pos+query_buffer.ptr-2)) && (pos != 0) &&
            (searchdata == bytestring(response_buffer.data[pos:(pos+query_buffer.ptr-2)]))
            return true
        end

        # Start searching
        # First the current response buffer
        match = backwards ?
                rsearch(bytestring(response_buffer.data[1:response_buffer.size]),searchdata,response_buffer.ptr - 1):
                response_buffer.ptr + 1 < response_buffer.size ?
                search(bytestring(response_buffer.data[1:response_buffer.size]),searchdata,response_buffer.ptr + 1): 0:-1

        if match != 0:-1
            seek(response_buffer,first(match)-1)
            return true
        end

        # Now search all the other buffers
        idx = hist.cur_idx
        found = false
        while true
            idx += backwards ? -1 : 1
            if !(0 < idx <= length(hist.history))
                break
            end
            match = backwards ? rsearch(hist.history[idx],searchdata):
                                search(hist.history[idx],searchdata);
            if match != 0:-1
                found = true
                truncate(response_buffer,0)
                write(response_buffer,hist.history[idx])
                seek(response_buffer,first(match)-1)
                break
            end
        end
        if found
            #if hist.cur_idx == length(hist.history)+1
            #    hist.last_buffer = copy(s.input_buffer)
            #end
            hist.cur_idx = idx
        end
        return found
    end

    function history_reset_state(hist::REPLHistoryProvider)
        hist.cur_idx = length(hist.history)+1
    end
    LineEdit.reset_state(hist::REPLHistoryProvider) = history_reset_state(hist)

    const julia_green = "\033[1m\033[32m"
    const color_normal = Base.color_normal

    function return_callback(repl,s)
        if position(s.input_buffer) != 0 && eof(s.input_buffer) &&
            (seek(s.input_buffer,position(s.input_buffer)-1);read(s.input_buffer,Uint8)=='\n')
            repl.consecutive_returns += 1
        else
            repl.consecutive_returns = 0
        end
        ast = parse_input_line(bytestring(copy(s.input_buffer.data)))
        if repl.consecutive_returns > 1 || !isa(ast,Expr) || (ast.head != :continue && ast.head != :incomplete)
            return true
        else
            return false
        end
    end

    function find_hist_file()
        if isfile(".julia_history2")
            return ".julia_history2"
        elseif haskey(ENV,"JULIA_HISTORY")
            return ENV["JULIA_HISTORY"]
        else
            @windows_only return ENV["AppData"]*"/julia/history2"
            @unix_only return ENV["HOME"]*"/.julia_history2"
        end
    end

    function send_to_backend(ast,req,rep)
        put!(req, (ast,1))
        (val, bt) = take!(rep)
    end

    have_color(s) = true

    function respond(f,d,main,req,rep)
        (s,buf,ok)->begin
            if !ok
                return transition(s,:abort)
            end
            line = takebuf_string(buf)
            if !isempty(line)
                reset(d)
                (val,bt) = send_to_backend(f(line),req,rep)
                if !ends_with_semicolon(line) || !is(bt,nothing)
                    print_response(d,val,bt,true,have_color(s))
                end
            end
            println(d.repl.t)
            reset_state(s)
            transition(s,main)
        end
    end

    import Terminals: raw!

    function reset(d::REPLDisplay)
        raw!(d.repl.t,false)
        print(Base.text_colors[:normal])
    end

    function setup_interface(d::REPLDisplay,req,rep;extra_repl_keymap=Dict{Any,Any}[])
        ###
        #
        # This function returns the main interface that describes the REPL
        # functionality, it is called internally by functions that setup a
        # Terminal-based REPL frontend, but if you want to customize your REPL
        # or embed the REPL in another interface, you may call this function
        # directly and append it to your interface.
        #
        # Usage:
        #
        # repl_channel,response_channel = RemoteRef(),RemoteRef()
        # start_repl_backend(repl_channel, response_channel)
        # setup_interface(REPLDisplay(t),repl_channel,response_channel)
        #
        ###

        ###
        # We setup the interface in two stages.
        # First, we set up all components (prompt,rsearch,shell,help)
        # Second, we create keymaps with appropriate transitions between them
        #   and assign them to the components
        #
        ###

        ############################### Stage I ################################

        repl = d.repl

        # This will provide completions for REPL and help mode
        replc = REPLCompletionProvider(repl)
        finalizer(replc,(replc)->close(f))

        # Set up the main Julia prompt
        main_prompt = Prompt("julia> ";
            # Copy colors from the prompt object
            prompt_color=repl.prompt_color,
            input_color=repl.input_color,
            keymap_func_data = repl,
            complete=replc,
            on_enter=s->return_callback(repl,s))

        main_prompt.on_done = respond(Base.parse_input_line,d,main_prompt,req,rep)

        # Setup help mode
        help_mode = Prompt("help> ",
            prompt_color = repl.help_color,
            input_color=repl.input_color,
            keymap_func_data = repl,
            complete = replc,
            # When we're done transform the entered line into a call to help("$line")
            on_done = respond(d,main_prompt,req,rep) do line
                parse("Base.@help $line")
            end)

        # Set up shell mode
        shell_mode = Prompt("shell> ";
            prompt_color = repl.shell_color,
            input_color=repl.input_color,
            keymap_func_data = repl,
            complete = ShellCompletionProvider(repl),
            # Transform "foo bar baz" into `foo bar baz` (shell quoting)
            # and pass into Base.repl_cmd for processing (handles `ls` and `cd`
            # special)
            on_done = respond(d,main_prompt,req,rep) do line
                Expr(:call, :(Base.repl_cmd), macroexpand(Expr(:macrocall,symbol("@cmd"),line)))
            end)

        ################################# Stage II #############################

        # Setup history
        # We will have a unified history for all REPL modes
        f = open(find_hist_file(),true,true,true,false,false)
        hp = hist_from_file(f,(Uint8=>Any)[uint8('\0') => main_prompt, uint8(';') => shell_mode, uint8('?') => help_mode, uint8('>') => main_prompt])
        history_reset_state(hp)
        main_prompt.hist = hp
        shell_mode.hist = hp
        help_mode.hist = hp

        (hkp,hkeymap) = LineEdit.setup_search_keymap(hp)

        # Canoniczlize user keymap input
        if isa(extra_repl_keymap,Dict)
            extra_repl_keymap = [extra_repl_keymap]
        end

        const repl_keymap = {
            ';' => function (s)
                if isempty(s) || position(LineEdit.buffer(s)) == 0
                    buf = copy(LineEdit.buffer(s))
                    transition(s,shell_mode)
                    LineEdit.state(s,shell_mode).input_buffer = buf
                    LineEdit.refresh_line(s)
                else
                    edit_insert(s,';')
                end
            end,
            '?' => function (s)
                if isempty(s) || position(LineEdit.buffer(s)) == 0
                    buf = copy(LineEdit.buffer(s))
                    transition(s,help_mode)
                    LineEdit.state(s,help_mode).input_buffer = buf
                    LineEdit.refresh_line(s)
                else
                    edit_insert(s,'?')
                end
            end,

            # Bracketed Paste Mode
            "\e[200~" => s->begin
                ps = LineEdit.state(s,LineEdit.mode(s))
                input = readuntil(ps.terminal,"\e[201~")[1:(end-6)]
                input = replace(input,'\r','\n')
                if position(LineEdit.buffer(s)) == 0
                    indent = Base.indentation(input)[1]
                    input = Base.unindent(input[(indent+1):end],indent)
                end
                buf = copy(LineEdit.buffer(s))
                edit_insert(buf,input)
                string = takebuf_string(buf)
                pos = 0
                while pos <= length(string)
                    oldpos = pos
                    ast, pos = Base.parse(string,pos; raise=false)
                    # Get the line and strip leading and trailing whitespace
                    line = strip(string[max(oldpos,1):min(pos-1,length(string))])
                    LineEdit.replace_line(s,line)
                    LineEdit.refresh_line(s)
                    if !isa(ast,Expr) || (ast.head != :continue && ast.head != :incomplete)
                        LineEdit.commit_line(s)
                        # This is slightly ugly but ok for now
                        terminal = LineEdit.terminal(s)
                        stop_reading(terminal)
                        raw!(terminal,false) && Terminals.Unix.disable_bracketed_paste(terminal)
                        LineEdit.mode(s).on_done(s,LineEdit.buffer(s),true)
                        raw!(terminal,true) && Terminals.Unix.enable_bracketed_paste(terminal)
                        start_reading(terminal)
                    else
                        break
                    end
                end
            end,
        }

        a = Dict{Any,Any}[hkeymap, repl_keymap, LineEdit.history_keymap(hp), LineEdit.default_keymap,LineEdit.escape_defaults]
        prepend!(a,extra_repl_keymap)
        @eval @LineEdit.keymap repl_keymap_func $(a)

        main_prompt.keymap_func = repl_keymap_func

        const mode_keymap = {
            '\b' => s->(isempty(s) ? transition(s,main_prompt) : LineEdit.edit_backspace(s) )
        }

        b = Dict{Any,Any}[hkeymap, mode_keymap, LineEdit.history_keymap(hp), LineEdit.default_keymap,LineEdit.escape_defaults]

        @eval @LineEdit.keymap mode_keymap_func $(b)

        shell_mode.keymap_func = help_mode.keymap_func = mode_keymap_func

        ModalInterface([main_prompt,shell_mode,help_mode,hkp])
    end

    run_frontend(repl::LineEditREPL,repl_channel,response_channel) = run_interface(repl.t,setup_interface(REPLDisplay(repl),repl_channel,response_channel))

    if isdefined(Base,:banner_color)
        banner(io,t) = banner(io,hascolor(t))
        banner(io,x::Bool) = print(io,x ? Base.banner_color : Base.banner_plain)
    else
        banner(io,t) = Base.banner(io)
    end

    function run_repl(t::TextTerminal)
        repl_channel = RemoteRef()
        response_channel = RemoteRef()
        start_repl_backend(repl_channel, response_channel)
        banner(t,t)
        run_frontend(LineEditREPL(t),repl_channel,response_channel)
    end

    type BasicREPL <: AbstractREPL
    end

    outstream(::BasicREPL) = STDOUT

    type StreamREPL <: AbstractREPL
        stream::IO
        prompt_color::String
        input_color::String
        answer_color::String
    end

    import Base.AsyncStream

    outstream(s::StreamREPL) = s.stream

    StreamREPL(stream::AsyncStream) = StreamREPL(stream,julia_green,Base.text_colors[:white],Base.answer_color())

    answer_color(r::LineEditREPL) = r.answer_color
    answer_color(r::StreamREPL) = r.answer_color
    answer_color(::BasicREPL) = Base.text_colors[:white]

    print_response(d::REPLDisplay,r::StreamREPL,args...) = print_response(d, r.stream,r, args...)
    print_response(d::REPLDisplay,r::LineEditREPL,args...) = print_response(d, r.t, r, args...)
    print_response(d::REPLDisplay,args...) = print_response(d,d.repl,args...)

    function run_repl(stream::AsyncStream)
    repl =
    @async begin
        repl_channel = RemoteRef()
        response_channel = RemoteRef()
        start_repl_backend(repl_channel,response_channel)
        StreamREPL_frontend(repl,repl_channel,response_channel)
    end
    repl
    end

    function ends_with_semicolon(line)
        match = rsearch(line,';')
        if match != 0
            for c in line[(match+1):end]
                if !isspace(c)
                    return c == '#'
                end
            end
            return true
        end
        return false
    end

    function run_frontend(repl::StreamREPL,repl_channel,response_channel)
        have_color = true
        banner(repl.stream,have_color)
        d = REPLDisplay(repl)
        while repl.stream.open
            if have_color
                print(repl.stream,repl.prompt_color)
            end
            print(repl.stream,"julia> ")
            if have_color
                print(repl.stream,repl.input_color)
            end
            line = readline(repl.stream)
            if !isempty(line)
                ast = Base.parse_input_line(line)
                if have_color
                    print(repl.stream,color_normal)
                end
                put!(repl_channel, (ast,1))
                (val, bt) = take!(response_channel)
                if !ends_with_semicolon(line)
                    print_response(d,val,bt,true,have_color)
                end
            end
        end
        # Terminate Backend
        put!(repl_channel,(nothing,-1))
    end

    function start_repl_server(port)
        listen(port) do server, status
            client = accept(server)
            run_repl(client)
        end
    end
end


text:
#!/usr/local/bin/julia
using Markdown

input = readall(STDIN)

cd(dirname(@__FILE__)) do
  open("log.txt", "w") do io
    println(io, input)
    println(io, "----")
  end
end

if input in ("\n", "")
  println("no input")
  exit(0)
end

try
  output = stringmime("text/html", Markdown.parse(input, flavor = :standard))

  println(output)

  cd(dirname(@__FILE__)) do
    open("log.txt", "a") do io
      println(io, output)
    end
  end

catch e
  open("log.txt", "a") do io
    showerror(io, e, catch_backtrace())
    println(io)
  end
#   rethrow()
end


text:
"Return a list of all indices in an n-dimensional hyperbolic cross."
function index_set_hyperbolic_cross(s, n, α = 1)
    I = _index_set_hyperbolic_cross(s, n, α)
    [tuple((1 .+ i)...) for i in I]
end

function _index_set_hyperbolic_cross(s, n, α = 1)
    if n == 1
        smax = floor(Int, s^(1/α))-1
        I = [[i] for i in 0:smax]
    else
        I = Array{Array{Int,1},1}()
        I_rec = _index_set_total_degree(s, n-1)
        for idx in I_rec
            for m in 0:floor(Int,s^(1/α)/prod(1+abs.(idx)))-1
                push!(I, [idx...; m])
            end
        end
        I
    end
end

"Return a list of all indices of total degree at most s, in n dimensions."
function index_set_total_degree(s, n)
    # We make a list of arrays because that is easier
    I = _index_set_total_degree(s, n)
    # and then we convert to tuples
    [tuple((1 .+ i)...) for i in I]
end

function _index_set_total_degree(s, n)
    if n == 1
        I = [[i] for i in 0:s]
    else
        I = Array{Array{Int,1},1}()
        I_rec = _index_set_total_degree(s, n-1)
        for idx in I_rec
            for m in 0:s-sum(abs.(idx))
                push!(I, [idx...; m])
            end
        end
        I
    end
end


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

using Markdown
using InteractiveUtils

# ╔═╡ 14f719cd-ffdc-4d9b-83f4-efd7e7b54c82
using Pkg; Pkg.activate(".")

# ╔═╡ c27879d8-9702-11eb-3050-f1c4477273fe
using PlotlyLight

# ╔═╡ e2dedee2-9702-11eb-1364-118bcaff2607
Plot(Config(x=1:50,y=randn(50), type="bar"))

# ╔═╡ Cell order:
# ╠═14f719cd-ffdc-4d9b-83f4-efd7e7b54c82
# ╠═c27879d8-9702-11eb-3050-f1c4477273fe
# ╠═e2dedee2-9702-11eb-1364-118bcaff2607


text:
@testset begin "dynamic selection"

    s = select(:x, :y => :z, :y => :w)

    # test Base.in
    @test :x in s
    @test (:x => :u) in s
    @test (:x => :u => :v) in s
    @test (:y => :z) in s
    @test (:y => :w) in s
    @test !(:y in s)
    @test !((:y => :r) in s)

    # test Base.getindex
    @test s[:x] == AllSelection()
    sub = s[:y]
    @test isa(sub, DynamicSelection)
    @test :z in sub
    @test :w in sub
    @test s[:u] == EmptySelection()

    # test set_subselection!
    set_subselection!(s, :y, select(:z))
    @test (:y => :z) in s
    @test !((:y => :w) in s)

    selection = select(:x)
    @test :x in selection
    subselection = select(:y)
    set_subselection!(selection, :x, subselection)
    @test (:x => :y) in selection
    @test !(:x in selection)
end

@testset begin "all selection"

    s = selectall()

    # test Base.in
    @test :x in s
    @test (:x => :u) in s

    # test Base.getindex
    @test s[:x] == AllSelection()
    @test s[:x => :y] == AllSelection()
end

max_stars_repo_name: dendisuhubdy/Gen

max_issues_repo_name: kumsh/Gen

max_forks_repo_name: kumsh/Gen


text:
test = """
    NNCB

    CH -> B
    HH -> N
    CB -> H
    NH -> C
    HB -> C
    HC -> B
    HN -> C
    NN -> C
    BH -> H
    NC -> B
    NB -> B
    BN -> B
    BB -> N
    BC -> B
    CC -> N
    CN -> C
    """

function pair_stats(s::AbstractString)::Dict
    pairs = [s[i:i+1] for i in 1:length(s) - 1]
    return Dict([p => count(==(p), pairs) for p in unique(pairs)])
end

function parse_input(io::IO)
    template = readline(io)
    readline(io)
    # rules0: between these -> insert this
    rules0 = split.(readlines(io), " -> ")
    # rules: this pair is created -> from these pairs
    rules = Dict{String, Vector{String}}()
    for (pat, ins) in rules0
        push!(get!(rules, string(pat[1], ins), []), pat)
        push!(get!(rules, string(ins, pat[2]), []), pat)
    end
    # first and last element never changes
    return (pair_stats(template), rules, last(template))
end

function step(pairs::AbstractDict, rules::AbstractDict)
    return Dict(a => sum([get(pairs, b, 0) for b in bs]) for (a, bs) in rules)
end

function run(polymer::AbstractDict, rules::AbstractDict, n::Int)::Dict
    for _ in 1:n
        polymer = step(polymer, rules)
    end
    return Dict(k => v for (k, v) in polymer if v > 0)
end

function element_count_range(polymer::AbstractDict, tail::AbstractChar)::Int
    counts = [[first(p) => c for (p, c) in polymer]; tail => 1]
    freqs = [sum([c for (elt, c) in counts if elt == unique_elt])
             for unique_elt in unique(first.(counts))]
    (min, max) = extrema(freqs)
    return max - min
end

q1(polymer, rules, tail) = element_count_range(run(polymer, rules, 10), tail)
q2(polymer, rules, tail) = element_count_range(run(polymer, rules, 40), tail)

let (polymer, rules, tail) = parse_input(IOBuffer(test))
    @assert pair_stats("NCNBCHB") == run(polymer, rules, 1)
    @assert pair_stats("NBCCNBBBCBHCB") == run(polymer, rules, 2)
    @assert pair_stats("NBBBCNCCNBBNBNBBCHBHHBCHB") == run(polymer, rules, 3)
    @assert q1(polymer, rules, tail) == 1588
    @assert q2(polymer, rules, tail) == 2188189693529
end

open("day14.in", "r") do io
    input = parse_input(io)
    println("Q1: ", q1(input...))
    println("Q2: ", q2(input...))
end

max_stars_repo_name: dhanak/competitive-coding

max_issues_repo_name: dhanak/competitive-coding

max_forks_repo_name: dhanak/competitive-coding


text:
@safetestset "isbase58" begin
    using Validator

    @test Validator.isbase58("BukQL") === true
    @test Validator.isbase58("3KMUV89zab") === true
    @test Validator.isbase58("91GHkLMNtyo98") === true
    @test Validator.isbase58("YyjKm3H") === true
    @test Validator.isbase58("Mkhss145TRFg") === true
    @test Validator.isbase58("7678765677") === true
    @test Validator.isbase58("abcodpq") === true
    @test Validator.isbase58("AAVHJKLPY") === true

    @test Validator.isbase58("0OPLJH") === false
    @test Validator.isbase58("IMKLP23") === false
    @test Validator.isbase58("KLMOmk986") === false
    @test Validator.isbase58("LL1l1985hG") === false
    @test Validator.isbase58("*MP9K") === false
    @test Validator.isbase58("Zm=8JBSWY3DP") === false
    @test Validator.isbase58(")()(=9292929MKL") === false
end


text:
function copy_user(X::Array{Int})
    x = Array{Int}(1:(size(X, 1)*size(X, 2)))
    for i in 1:(size(X, 1)*size(X, 2))
        x[i] = X[i];
    end
    return reshape(x, (size(X, 2), size(X, 1)))
end

print(copy_user(Vector(2:5)))

print(copy_user([[2, 7, 8] [6, 3, 4] [9, 1, 0]]))


text:
using HDF5
import Base: findall

function get_fnames_and_types(filepath)
    filenames = first.(splitext.(last.(split.(filepath,'/'))))
    fileext = last.(splitext.(last.(split.(filepath,'/'))))

    return filenames,fileext
end
"""
    ec_grab(measurement::AbstractString; dir::AbstractString = "./Data/EC")

Fetch the CV and OCP Measurements from your dir.

eg.
    julia> v,c = ec_grab("001-001",dir = "..")\n\t

You can also grab the DataFrame from ec_list(). Be sure that the experiments in ec_list() all are either CV or OCP or another experiment. Otherwise this will result in an error.
    julia>v=ec_grab(ec_list(measurement="OCP"))

    julia>v=ec_grab(ec_list(measurement="CV"))
"""
function ec_grab(measurement::AbstractString; dir::AbstractString = "./Data/EC")
    if occursin("EC-",measurement) == false
        _measurement = "EC-"*measurement
    else
        _measurement = measurement
    end

        dirs_files = readdir(dir,join=true)
        files = [dirs_files[i] for i in 1:length(dirs_files) if isdir(dirs_files[i]) == false]
        filenames,fileext =get_fnames_and_types(files)

        idxs = findall(x-> x==(_measurement),filenames)

        if isempty(idxs) == true
            error("There is no Measurement named $(_measurement) in $(dir).")
        else
            idx_file = idxs[1]

            if fileext[idx_file] == ".h5"
                _measurement =  get_measurement(files[idx_file])

                if _measurement == "OCP"
                    _volt = h5open(files[idx_file]) do fid
                                read(fid["OCP"]["Data"])
                            end

                    return _volt

                elseif _measurement == "CV"
                    data = h5open(files[idx_file]) do fid
                        read(fid["CV"]["Data"])
                        end
                    _volt = data[:,1]
                    _curr = data[:,2]

                    return _volt,_curr
                    end

            elseif fileext[idx_file] == ".csv"
                data = readdlm(files[idx_file])
                _volt = data[:,1]
                _curr = data[:,2]
                    return _volt,_curr

            elseif fileext[idx_file] == ".txt"
                data = readdlm(files[idx_file])
                if typeof(data[1]) == String
                    _volt = data[:,1]
                    _curr = data[:,2]
                    return _volt,_curr
                else
                    _volt = data[2:end,1]
                    _curr = data[2:end,2]
                    return _volt,_curr
                end
            end
        end

end

function ec_grab(df::DataFrame)

    path = df[:,end]
    fileext= get_fnames_and_types(path)[2]

    data=[]
    for i in 1:first(size(df))
        if fileext[i] == ".h5"
            _measurement = get_measurement(path[i])
            _data = h5open(path[i]) do fid
                read(fid[_measurement]["Data"])
            end
            push!(data,_data)
        elseif fileext[i] == ".csv"
            _data = readdlm(path[i])
            push!(data,_data)
        elseif fileext[i] == ".txt"
            _data = readdlm(path[i])[:,1:2]
            push!(data,_data)
        end
    end

    if all(last.(size.(data)).== last(size(data[1])))
        if (last(size(data[1])) == 1) || (try size(data[1])[2] catch end === nothing)
            _dim1 = [data[i][:,1] for i in 1:length(data)]
            @show data[1]
            return _dim1
        elseif last(size(data[1])) == 2
            _dim1 = [data[i][:,1] for i in 1:length(data)]
            _dim2 = [data[i][:,2] for i in 1:length(data)]
            return _dim1,_dim2
        end
    else
        error("The Files you selected in ec_list() are not compatible.")
    end

end

"""
Check if measurement(.h5 file) is a CV or OCP or else.
    get_measurement(file::AbstractString)

e.g.
    julia> get_measurement("./EC-CV.h5")
    "CV"
"""
function get_measurement(file::AbstractString)

    _measurement =  h5open(file) do fid
        keys(fid)[1]
    end
end

"""
Check if current range from CV or CM measurment(.h5 file).
    get_curr_range(file::AbstractString)

e.g.
    julia> get_curr_range("./EC-CV.h5")
    1e-6
"""
function get_curr_range(file::AbstractString)
    _measurement = get_measurement(file)

    _curr_range =  h5open(file) do fid
        HDF5.attributes(fid[_measurement]["Data"])["Current Range [A]"] |> read
    end
end

"""
Get Parameters from h5.file
    get_params(file::AbstractString)

e.g.
    julia> get_params("./EC-CV.h5")
    5×2 Matrix{Any}:
     "Sample Rate [Hz]"        10.0
     "Scan Rate [mV/s]"        50.0
     "Current Range [A]"        1.0e-6
     "Current Filter [Hz]"     10.0
     "Bandwidth Filter [Hz]"  100.0
"""
function get_params(file::AbstractString)
    _measurement =get_measurement(file)
    if _measurement == "CV"

        attributes = [
            "Sample Rate [Hz]",
            "Scan Rate [mV/s]",
            "Current Range [A]",
            "Current Filter [Hz]",
            "Bandwidth Filter [Hz]"
        ]

        values = [HDF5.attributes(h5open(file)["CV"]["Data"])[attributes[i]] |> read for i in 1:length(attributes)]
        close(h5open(file))
        return [attributes values]

    elseif _measurement == "OCP"

        attribute = "Sample Rate [Hz]"

        value = HDF5.attributes(h5open(file)["OCP"]["Data"])[attribute] |> read
        close(h5open(file))

        return [attribute value]

    end

end

"""
Returns the List of all Measurements in dir. Use kwargs to filter the df.
Use show_na to show older measurements without parameters.\t

    function ec_list(;   dir::AbstractString = "./test/Data/EC",
        inexact::AbstractString  = "",
        measurement::AbstractString = "",
        scan_rate = "",
        show_na=false
   )

e.g.\t\n
    julia> ec_list()
    2×4 DataFrame
    Row │ Filename     Measurement  ScanRate  Path
        │ String       String       Any       String
   ─────┼─────────────────────────────────────────────────────────────────── \n
      1 │ EC-test_CV   CV           50.0      ./test/Data/EC/EC-test_CV.h5
      2 │ EC-test_OCP  OCP          n/A       ./test/Data/EC/EC-test_OCP.h5

    julia> ec_list(show_na=true)
    4×4 DataFrame
    Row │ Filename     Measurement  ScanRate  Path
        │ String       String       Any       String
   ─────┼───────────────────────────────────────────────────────────────────\n
      1 │ EC-001-001   n/A          n/A       ./test/Data/EC/EC-001-001.txt
      2 │ EC-XXX-XXX   n/A          n/A       ./test/Data/EC/EC-XXX-XXX.csv
      3 │ EC-test_CV   CV           50.0      ./test/Data/EC/EC-test_CV.h5
      4 │ EC-test_OCP  OCP          n/A       ./test/Data/EC/EC-test_OCP.h5

    julia> ec_list(inexact="001-001",show_na=true)
    1×4 DataFrame
    Row │ Filename    Measurement  ScanRate  Path
        │ String      String       Any       String
   ─────┼──────────────────────────────────────────────────────────────────\n
      1 │ EC-001-001  n/A          n/A       ./test/Data/EC/EC-001-001.txt

    julia> ec_list(scan_rate=50)
    1×4 DataFrame
    Row │ Filename    Measurement  ScanRate  Path
        │ String      String       Any       String
   ─────┼─────────────────────────────────────────────────────────────────\n
      1 │ EC-test_CV  CV           50.0      ./test/Data/EC/EC-test_CV.h5
"""
function ec_list(;   dir::AbstractString = "./test/Data/EC",
                     inexact::AbstractString  = "",
                     measurement::AbstractString = "",
                     scan_rate = missing ,
                     #pump::AbstractString = "",
                     show_na = false
                    )

    dirs_files = readdir(dir,join=true)
    files = [dirs_files[i] for i in 1:length(dirs_files) if isdir(dirs_files[i]) == false]
    if Sys.isapple() == true
        filenames = first.(splitext.(last.(split.(files,'/'))))
        filetype = last.(splitext.(last.(split.(files,'/'))))
    elseif Sys.iswindows() == true
        filenames = first.(splitext.(last.(split.(files,"\\"))))
        filetype = last.(splitext.(last.(split.(files,"\\"))))
    end

    measurements = [filetype[i] == ".h5" ? get_measurement(files[i]) : "n/A" for i in 1:length(files)]
    scan_rates = [(filetype[i] == ".h5" && get_measurement(files[i]) == "CV") ? get_params(files[i])[2,2] : "n/A" for i in 1:length(files)]
    _df = DataFrame(Filename = filenames, Measurement = measurements, ScanRate= scan_rates, Path=files)

   if show_na == false
        _df=_df[(_df.Measurement .!== "n/A"),:]
    if inexact !== ""
        _df=_df[(occursin.(inexact,_df.Filename) .== true ),:]
    end
    if measurement !== ""
        _df = _df[(_df.Measurement .== measurement),:]
    end
    if scan_rate !== missing
        _df = _df[(_df.ScanRate .== scan_rate),:]
    end
    return _df

   elseif show_na == true
    if inexact !== ""
        _df=_df[(occursin.(inexact,_df.Filename) .== true ),:]
    end
    if measurement !== ""
        _df = _df[((_df.Measurement .== measurement) .| (_df.Measurement .== "n/A")),:]
    end
    if scan_rate !== missing
        _df = _df[(((_df.ScanRate .== scan_rate) .| (_df.ScanRate .== "n/A")) .& (_df.Measurement .== "CV" )),:]
    end
        return _df

    end

end

max_stars_repo_name: lambcount/EC

max_issues_repo_name: lambcount/EC

max_forks_repo_name: lambcount/EC


text:
using CUSPARSE
using CUDArt
using Base.Test

m = 25
n = 35
k = 10
blockdim = 5

@testset "doti" begin
    @testset for elty in [Float32,Float64,Complex64,Complex128]
        x = sparsevec(rand(1:m,k), rand(elty,k), m)
        y = rand(elty,m)
        d_x = CudaSparseVector(x)
        d_y = CudaArray(y)
        ddot = CUSPARSE.doti(d_x,d_y,'O')
        #compare
        dot = zero(elty)
        for i in 1:length(x.nzval)
            dot += x.nzval[i] * y[x.nzind[i]]
        end
        @test ddot ≈ dot
    end
end

@testset "dotci" begin
    @testset for elty in [Complex64,Complex128]
        x = sparsevec(rand(1:m,k), rand(elty,k), m)
        y = rand(elty,m)
        d_x = CudaSparseVector(x)
        d_y = CudaArray(y)
        ddot = CUSPARSE.dotci(d_x,d_y,'O')
        #compare
        dot = zero(elty)
        for i in 1:length(x.nzval)
            dot += conj(x.nzval[i]) * y[x.nzind[i]]
        end
        @test ddot ≈ dot
    end
end


text:
module Euler
  export EulerLaw, γ

  import ..Atum
  using ..Atum: avg, logavg, roe_avg
  using StaticArrays
  using LinearAlgebra: I

  struct EulerLaw{γ, FT, D, S} <: Atum.AbstractBalanceLaw{FT, D, S}
    function EulerLaw{FT, D}(; γ = 7 // 5) where {FT, D}
      S = 2 + D
      new{FT(γ), FT, D, S}()
    end
  end

  γ(::EulerLaw{_γ}) where {_γ} = _γ

  function varsindices(law::EulerLaw)
    S = Atum.numberofstates(law)
    ix_ρ = 1
    ix_ρu⃗ = StaticArrays.SUnitRange(2, S - 1)
    ix_ρe = S
    return ix_ρ, ix_ρu⃗, ix_ρe
  end

  function unpackstate(law::EulerLaw, q)
    ix_ρ, ix_ρu⃗, ix_ρe = varsindices(law)
    @inbounds q[ix_ρ], q[ix_ρu⃗], q[ix_ρe]
  end

  function pressure(law::EulerLaw, ρ, ρu⃗, ρe)
    (γ(law) - 1) * (ρe - ρu⃗' * ρu⃗ / 2ρ)
  end
  function energy(law::EulerLaw, ρ, ρu⃗, p)
    p / (γ(law) - 1) + ρu⃗' * ρu⃗ / 2ρ
  end
  function soundspeed(law::EulerLaw, ρ, p)
    sqrt(γ(law) * p / ρ)
  end
  function soundspeed(law::EulerLaw, ρ, ρu⃗, ρe)
    soundspeed(law, ρ, pressure(law, ρ, ρu⃗, ρe))
  end

  function Atum.flux(law::EulerLaw, q, _)
    ρ, ρu⃗, ρe = unpackstate(law, q)

    u⃗ = ρu⃗ / ρ
    p = pressure(law, ρ, ρu⃗, ρe)

    fρ = ρu⃗
    fρu⃗ = ρu⃗ * u⃗' + p * I
    fρe = u⃗ * (ρe + p)

    hcat(fρ, fρu⃗, fρe)
  end

  function Atum.wavespeed(law::EulerLaw, n⃗, q, x⃗)
    ρ, ρu⃗, ρe = unpackstate(law, q)

    u⃗ = ρu⃗ / ρ
    abs(n⃗' * u⃗) + soundspeed(law, ρ, ρu⃗, ρe)
  end

  function Atum.surfaceflux(::Atum.RoeFlux, law::EulerLaw, n⃗, x⃗, q⁻, q⁺)
    f⁻ = Atum.flux(law, q⁻, x⃗)
    f⁺ = Atum.flux(law, q⁺, x⃗)

    ρ⁻, ρu⃗⁻, ρe⁻ = unpackstate(law, q⁻)
    u⃗⁻ = ρu⃗⁻ / ρ⁻
    e⁻ = ρe⁻ / ρ⁻
    p⁻ = pressure(law, ρ⁻, ρu⃗⁻, ρe⁻)
    h⁻ = e⁻ + p⁻ / ρ⁻
    c⁻ = soundspeed(law, ρ⁻, p⁻)

    ρ⁺, ρu⃗⁺, ρe⁺ = unpackstate(law, q⁺)
    u⃗⁺ = ρu⃗⁺ / ρ⁺
    e⁺ = ρe⁺ / ρ⁺
    p⁺ = pressure(law, ρ⁺, ρu⃗⁺, ρe⁺)
    h⁺ = e⁺ + p⁺ / ρ⁺
    c⁺ = soundspeed(law, ρ⁺, p⁺)

    ρ = sqrt(ρ⁻ * ρ⁺)
    u⃗ = roe_avg(ρ⁻, ρ⁺, u⃗⁻, u⃗⁺)
    h = roe_avg(ρ⁻, ρ⁺, h⁻, h⁺)
    c = roe_avg(ρ⁻, ρ⁺, c⁻, c⁺)

    uₙ = u⃗' * n⃗

    Δρ = ρ⁺ - ρ⁻
    Δp = p⁺ - p⁻
    Δu⃗ = u⃗⁺ - u⃗⁻
    Δuₙ = Δu⃗' * n⃗

    c⁻² = 1 / c^2
    w1 = abs(uₙ - c) * (Δp - ρ * c * Δuₙ) * c⁻² / 2
    w2 = abs(uₙ + c) * (Δp + ρ * c * Δuₙ) * c⁻² / 2
    w3 = abs(uₙ) * (Δρ - Δp * c⁻²)
    w4 = abs(uₙ) * ρ

    fp_ρ = (w1 + w2 + w3) / 2
    fp_ρu⃗ = (w1 * (u⃗ - c * n⃗) +
             w2 * (u⃗ + c * n⃗) +
             w3 * u⃗ +
             w4 * (Δu⃗ - Δuₙ * n⃗)) / 2
    fp_ρe = (w1 * (h - c * uₙ) +
             w2 * (h + c * uₙ) +
             w3 * u⃗' * u⃗ / 2 +
             w4 * (u⃗' * Δu⃗ - uₙ * Δuₙ)) / 2

    (f⁻ + f⁺)' * n⃗ / 2 - SVector(fp_ρ, fp_ρu⃗..., fp_ρe)
  end

  function Atum.twopointflux(::Atum.EntropyConservativeFlux,
                             law::EulerLaw,
                             q₁, _, q₂, _)
      FT = eltype(law)
      ρ₁, ρu⃗₁, ρe₁ = unpackstate(law, q₁)
      ρ₂, ρu⃗₂, ρe₂ = unpackstate(law, q₂)

      u⃗₁ = ρu⃗₁ / ρ₁
      p₁ = pressure(law, ρ₁, ρu⃗₁, ρe₁)
      b₁ = ρ₁ / 2p₁

      u⃗₂ = ρu⃗₂ / ρ₂
      p₂ = pressure(law, ρ₂, ρu⃗₂, ρe₂)
      b₂ = ρ₂ / 2p₂

      ρ_avg = avg(ρ₁, ρ₂)
      u⃗_avg = avg(u⃗₁, u⃗₂)
      b_avg = avg(b₁, b₂)

      u²_avg = avg(u⃗₁' * u⃗₁, u⃗₂' * u⃗₂)
      ρ_log = logavg(ρ₁, ρ₂)
      b_log = logavg(b₁, b₂)

      fρ = u⃗_avg * ρ_log
      fρu⃗ = u⃗_avg * fρ' + ρ_avg / 2b_avg * I
      fρe = (1 / (2 * (γ(law) - 1) * b_log) - u²_avg / 2) * fρ + fρu⃗ * u⃗_avg

      hcat(fρ, fρu⃗, fρe)
  end
end


text:
"""
    HTTPCompatAdapter <: ProtocolAdapter

An adapter that works for HTTP API v1.x.
"""
mutable struct HTTPCompatAdapter <: HTTPAdapterBase
    server::String
    output_channel::OutputChannel
    sessionKey::String
    closed::Bool
end

HTTPCompatAdapter() = HTTPCompatAdapter("", OutputChannel(DEFAULT_BUFFER_SIZE), "", false)

function loop(adp::HTTPCompatAdapter, server, qq, authKey; poll_interval = 1, fetch_count = 8)
    adp.server = server
    # Verify
    resp = post_restful("http://$server/auth", (), (; authKey))
    adp.sessionKey = resp[:session]
    # Bind
    resp = post_restful("http://$server/verify", (), (; adp.sessionKey, qq))
    resp[:code] == 0 || throw(AdapterConectionFailed(resp[:code], resp[:msg]))

    put!(adp.output_channel, resp)
    loop_body(adp, qq, poll_interval, fetch_count, Commands.response_type_compat, SESSION_KEY_IN_BODY)
end

Commands.response_type_compat(msg) = Commands.response_type(msg)

function send(adp::HTTPCompatAdapter, cmd::Commands.AbstractCommand; kwarg...)
    send(adp, cmd, Commands.response_type_compat; session_key_position = SESSION_KEY_IN_BODY, kwarg...)
end

is_adapter_compatibile(::HTTPCompatAdapter, version) = version.major == 1


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

using Markdown
using InteractiveUtils

# ╔═╡ eb3b3be0-67d1-11eb-144b-2f9581950503
md"""
#### Math can be easy!
#### Just write simple steps.
"""

# ╔═╡ 66294030-67d3-11eb-25b1-8d25703546db
md"""
#### sqrt
##### ??? What is this?
####  Pluto 🛡 to the rescue!
#####  Doing Things can be objects.
"""

# ╔═╡ f4835c0a-67db-11eb-134c-85092b5a1ced
md"""
##### Looking for a short
##### cut across a square 🔲?
##### We call it the root 🥕.
"""

# ╔═╡ 1976e31c-67df-11eb-2d4a-1d8e612a99a9

md"""
##### sqrt() is a Doing Thing
##### called a function.
##### It can tells how far it
##### is across a 1 by 1 🔲.
"""

# ╔═╡ 84554b8a-67ce-11eb-2e10-7550cbde4485
md"""
##### Hard to think of sqrt?
##### I bet you will never forget
## 🔲🥕
##### Lets create a new Doing Thing.
"""

# ╔═╡ f8c8c828-67d0-11eb-16b6-f384a4a4a5df
🔲🥕 = sqrt

# ╔═╡ 279bdd4e-67e2-11eb-1adb-b98afc7e812e
md"""

##### Yes, it works!
##### Will you remember?
### 🔲🥕  $~~~~or~~~~~$   sqrt

"""

# ╔═╡ 1c3a7e52-67d1-11eb-0c5a-a5b36b4fc989
🔲🥕(2)

# ╔═╡ Cell order:
# ╠═eb3b3be0-67d1-11eb-144b-2f9581950503
# ╠═66294030-67d3-11eb-25b1-8d25703546db
# ╠═f4835c0a-67db-11eb-134c-85092b5a1ced
# ╠═1976e31c-67df-11eb-2d4a-1d8e612a99a9
# ╠═84554b8a-67ce-11eb-2e10-7550cbde4485
# ╠═f8c8c828-67d0-11eb-16b6-f384a4a4a5df
# ╠═279bdd4e-67e2-11eb-1adb-b98afc7e812e
# ╠═1c3a7e52-67d1-11eb-0c5a-a5b36b4fc989

max_stars_repo_name: paradocs/my-first-binder

max_issues_repo_name: paradocs/my-first-binder

max_forks_repo_name: paradocs/my-first-binder


text:
# ---
# title: 1223. Dice Roll Simulation
# id: problem1223
# author: Tian Jun
# date: 2020-10-31
# difficulty: Medium
# categories: Dynamic Programming
# link: <https://leetcode.com/problems/dice-roll-simulation/description/>
# hidden: true
# ---
#
# A die simulator generates a random number from 1 to 6 for each roll. You
# introduced a constraint to the generator such that it cannot roll the number
# `i` more than `rollMax[i]` (1-indexed) **consecutive** times.
#
# Given an array of integers `rollMax` and an integer `n`, return the number of
# distinct sequences that can be obtained with exact `n` rolls.
#
# Two sequences are considered different if at least one element differs from
# each other. Since the answer may be too large, return it modulo `10^9 + 7`.
#
#
#
# **Example 1:**
#
#
#
#     Input: n = 2, rollMax = [1,1,2,2,2,3]
#     Output: 34
#     Explanation: There will be 2 rolls of die, if there are no constraints on the die, there are 6 * 6 = 36 possible combinations. In this case, looking at rollMax array, the numbers 1 and 2 appear at most once consecutively, therefore sequences (1,1) and (2,2) cannot occur, so the final answer is 36-2 = 34.
#
#
# **Example 2:**
#
#
#
#     Input: n = 2, rollMax = [1,1,1,1,1,1]
#     Output: 30
#
#
# **Example 3:**
#
#
#
#     Input: n = 3, rollMax = [1,1,1,2,2,3]
#     Output: 181
#
#
#
#
# **Constraints:**
#
#   * `1 <= n <= 5000`
#   * `rollMax.length == 6`
#   * `1 <= rollMax[i] <= 15`
#
#
## @lc code=start
using LeetCode

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


text:
@testset "Spectrum" begin
    wave = range(1e4, 5e4, length = 1000)
    sigma = randn(size(wave))
    sigma[7] = 1
    sigma[134] = 0.1
    flux = 100 .± sigma
    flux[7] = 1000 ± 1
    flux[134] = 1 ± 0.1

    spec = spectrum(wave, flux, name = "test spectrum")

    @test spec.wave == wave
    @test size(spec) === (1000,)
    @test length(spec) == 1000
    @test maximum(spec) == 1000 ± 1
    @test minimum(spec) == 1 ± 0.1
    @test argmax(spec) == 7
    @test argmin(spec) == 134
    @test findmax(spec) == (1000 ± 1, 7)
    @test findmin(spec) == (1 ± 0.1, 134)
    @test spec.flux == flux
    @test Measurements.uncertainty.(spec.flux) ≈ sigma

    flux_trimmed = flux[200:800]
    @test_throws AssertionError spectrum(wave, flux_trimmed)
    expected = """
    Spectrum(Float64, Measurement{Float64})
      name: test spectrum"""
    @test sprint(show, spec) == expected
    @test spec.name == "test spectrum"
end

@testset "Unitful Spectrum" begin
    wave = range(1e4, 5e4, length = 1000)
    sigma = randn(size(wave))
    sigma[7] = 1
    sigma[134] = 0.1
    flux = 100 .± sigma
    flux[7] = 1000 ± 1
    flux[134] = 1 ± 0.1

    funit = u"W/m^2/angstrom"
    spec = spectrum(wave * u"angstrom", flux * funit, name = "test")

    @test spec.wave ≈ wave * u"angstrom"

    @test size(spec) === (1000,)
    @test length(spec) == 1000
    @test maximum(spec) == (1000 ± 1) * funit
    @test minimum(spec) == (1 ± 0.1) * funit
    @test argmax(spec) == 7
    @test argmin(spec) == 134
    @test findmax(spec) == ((1000 ± 1) * funit, 7)
    @test findmin(spec) == ((1 ± 0.1) * funit, 134)
    @test spec.name == "test"

    # Test stripping
    w_unit, f_unit = unit(spec)
    @test w_unit == u"angstrom"
    @test f_unit == u"W/m^2/angstrom"

    strip_spec = ustrip(spec)
    @test strip_spec.wave == ustrip.(spec.wave)
    @test strip_spec.flux == ustrip.(spec.flux)
    @test strip_spec.meta == spec.meta
    expected = """
    Spectrum(Quantity{Float64,𝐋,Unitful.FreeUnits{(Å,),𝐋,nothing}}, Quantity{Measurement{Float64},𝐌*𝐋^-1*𝐓^-3,Unitful.FreeUnits{(Å^-1, m^-2, W),𝐌*𝐋^-1*𝐓^-3,nothing}})
      name: test"""
    @test sprint(show, spec) == expected
end

@testset "Arithmetic" begin
    wave = range(1e4, 5e4, length = 1000)
    sigma = randn(size(wave))
    flux = 100 .± sigma

    spec = spectrum(wave, flux, name = "test spectrum")

    # Scalars/ vectors
    values = [10, randn(size(spec))]
    for A in values
        # addition
        s = spec + A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .+ A

        # subtraction
        s = spec - A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .- A

        # multiplication
        s = spec * A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .* A

        # division
        s = spec / A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux ./ A
    end

    # Other spectra
    spec2 = deepcopy(spec)

    s = spec + spec2
    @test s.wave == spec.wave
    @test s.flux ≈ 2 .* spec.flux

    s = spec - spec2
    @test s.wave == spec.wave
    @test s.flux ≈ zeros(size(spec))

    s = spec * spec2
    @test s.wave == spec.wave
    @test s.flux ≈ spec.flux.^2

    s = spec / spec2
    @test s.wave == spec.wave
    @test s.flux ≈ ones(size(spec))

    spec = spectrum(spec.wave * u"cm", spec.flux * u"W/m^2/cm", name = "test unitfulspectrum")

    # Scalars/ vectors
    for A in [10u"W/m^2/cm", randn(size(spec))u"W/m^2/cm"]
        # addition
        s = spec + A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .+ A

        # subtraction
        s = spec - A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .- A
    end

    for A in [10, randn(size(spec))]
        # multiplication
        s = spec * A
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux .* A

        # division
        s = spec / 10
        @test s.wave == spec.wave
        @test s.flux ≈ spec.flux ./ 10
    end
end


text:
# Define a function
f(x,y) = exp(-(x-y)^2)

# Evaluate the function
f(0.2, 1) # Returns 0.5272924240430484

# Evaluate with a Float32
f(0.2f0, 1) #Returns 0.5272924f0

@code_typed f(0.2, 1)

f(im*0.2, 1f0) # 0.3526677021528701 + 0.1491055129382033im

# Infer which function is called
@which exp(0.2f0)
@which exp(im*0.2f0)

# Simple definition of a data structure
using Base: Number
struct MyNumber <: Number end

# Dispatch on the definition
Base.exp(t::MyNumber) = println("Foo")

# Execute
a = MyNumber()
exp(a) # Returns "Foo"

# Import Linear Algebra
using LinearAlgebra

# Add a function dispatch on two vectors
f(x::AbstractVector, y::AbstractVector) = exp.(-(x-y).^2)

# Define two random vectors
x = randn(10); y = randn(10);

f(x,y) # Returns a vector

methods(f) # Show the available methods corresponding to f


text:
module blueqat

    include("./backends/backend_functions.jl")
    using Reexport
    @reexport using .BackendFunctions

    export QuantumCircuit, Input, apply!, get_device, execute

    devices = Dict([("UndirectedGraph", UndirectedGraphModel),
                    ("StateVector", StateVectorModel)
                        ])

    mutable struct QuantumCircuit
        _n_qreg::Int64
        _n_creg::Int64
        _gates::Array{Gate,1}
    end

    mutable struct Input
        _n_qreg::Int64
        _init_state::Array{Complex, 1}
        _gates::Array{Gate,1}
    end

    function get_device(name::String, x...)
        return devices[name](x)
    end

    function QuantumCircuit(n_qubit::Int64, n_creg::Int64)
        _n_qubit = n_qubit
        _n_creg = n_creg
        return QuantumCircuit(_n_qubit, _n_creg, [])
    end

    function Input(init_state::Array{Complex, 1})
        _init_state = init_state
        _n_qreg = length(init_state)
        return Input(n_qreg, init_state, [])
    end

    function Input(n_qreg::Int)
        _init_state = zeros(Complex, n_qreg)
        _n_qreg = n_qreg
        return Input(n_qreg, _init_state, [])
    end

    function apply!(circuit::QuantumCircuit, gate::Gate)
        push!(circuit._gates, gate)
        return circuit
    end

    function apply!(circuit::QuantumCircuit, gates::Array{T,1}) where T <: Gate
        for gate in gates
            push!(circuit._gates, gate)
        end
        return circuit
    end

    function apply!(circuit::QuantumCircuit, gates::Array{Any,1})
        for gate in gates
            if typeof(gate) <: Gate
                push!(circuit._gates, gate)
            else
                error("Incorrect gates.")
            end
        end
        return circuit
    end

    function apply!(circuit::Input, gate::Gate)
        push!(circuit._gates, gate)
        return circuit
    end

    function apply!(circuit::Input, gates::Array{T,1}) where T <: Gate
        for gate in gates
            push!(circuit._gates, gate)
        end
        return circuit
    end

    function apply!(circuit::Input, gates::Array{Any,1})
        for gate in gates
            if typeof(gate) <: Gate
                push!(circuit._gates, gate)
            else
                error("Incorrect gates.")
            end
        end
        return circuit
    end

    function execute(circuit::QuantumCircuit, device::T, args...) where T <: Device
        result = execute_backend(circuit._n_qreg::Int, circuit._gates::Array{Gate,1}, device, args)
        return result
    end

    function execute(circuit::Input, device::T, args...) where T <: Device
        result = execute_backend(circuit._n_qreg::Int, circuit._gates::Array{Gate,1}, device, args)
        return result
    end

end

max_stars_repo_name: terasakisatoshi/blueqat-julia

max_issues_repo_name: terasakisatoshi/blueqat-julia

max_forks_repo_name: terasakisatoshi/blueqat-julia


text:
function splintwp!(n, wp, f, g)
    # integer, intent(in) :: n
    # real(8), intent(in) :: wp[4,n),f(n)
    # real(8), intent(out) :: g(n)
    g[1] = 0.0
    g[2] = wp[1,2]*f[1] + wp[2,2]*f[2] + wp[3,2]*f[3] + wp[4,2]*f[4]
    for i in 3:n-1
      g[i] = g[i-1] + wp[1,i]*f[i-2] + wp[2,i]*f[i-1] + wp[3,i]*f[i] + wp[4,i]*f[i+1]
    end
    g[n] = g[n-1] + wp[1,n]*f[n-3] + wp[2,n]*f[n-2] + wp[3,n]*f[n-1] + wp[4,n]*f[n]
    return
end


text:
a::Int32 = a
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
ccall(:memset, Ptr{Void}, (Ptr{Void}, Int32, Csize_t),
    px, a, n*selty)


text:
import DataFrames.skipmissing

include("clusters.jl")
include("planetary_catalog.jl")
include("optimization.jl")

sim_param = setup_sim_param_model()

##### To start saving the model iterations in the optimization into a file:

add_param_fixed(sim_param,"num_targets_sim_pass_one",78005)

model_name = "Clustered_P_R_broken_R_optimization_blocks"
optimization_number = "_random"*ARGS[1] #if want to run on the cluster with random initial active parameters: "_random"*ARGS[1]
use_KS_or_AD = "KS" #'KS' or 'AD'
AD_mod = true
Kep_or_Sim = "Sim" #'Kep' or 'Sim'
max_evals = 300
num_evals_weights = 20
dists_exclude = [7,8,13] #Int64[] if want to include all distances
cycles = 3
Pop_per_param = 4

file_name = model_name*optimization_number*"_targs"*string(get_int(sim_param,"num_targets_sim_pass_one"))*"_evals"*string(max_evals)*"_cycles"*string(cycles)*".txt"
f = open(file_name, "w")
println(f, "# All initial parameters:")
write_model_params(f, sim_param)

##### To run the same model multiple times to see how it compares to a simulated catalog with the same parameters:

srand(1234) #to have the same reference catalog and simulated catalogs for calculating the weights

#To generate a simulated catalog to fit to:
cat_phys = generate_kepler_physical_catalog(sim_param)
cat_phys_cut = ExoplanetsSysSim.generate_obs_targets(cat_phys,sim_param)
cat_obs = observe_kepler_targets_single_obs(cat_phys_cut,sim_param)
summary_stat_ref = calc_summary_stats_model(cat_obs,sim_param)

#To simulate more observed planets for the subsequent model generations:
add_param_fixed(sim_param,"max_incl_sys",80.0) #degrees; 0 (deg) for isotropic system inclinations; set closer to 90 (deg) for more transiting systems

active_param_true, weights, target_fitness, target_fitness_std = compute_weights_target_fitness_std_perfect_model(num_evals_weights, use_KS_or_AD ; AD_mod=AD_mod, weight=true, dists_exclude=dists_exclude, save_dist=true)

##### To set up a block-optimization routine:

#To specify the blocks of random parameters, their search ranges, and to draw initial values:
active_param_keys = [["log_rate_clusters", "log_rate_planets_per_cluster", "num_mutual_hill_radii", "power_law_P", "sigma_logperiod_per_pl_in_cluster"], ["break_radius", "mr_power_index", "power_law_r1", "power_law_r2", "sigma_log_radius_in_cluster"], ["sigma_hk", "sigma_incl", "sigma_incl_near_mmr"]]
    #[["log_rate_clusters", "log_rate_planets_per_cluster", "num_mutual_hill_radii", "power_law_P", "sigma_logperiod_per_pl_in_cluster"], ["break_radius", "mr_power_index", "power_law_r1", "power_law_r2", "sigma_log_radius_in_cluster"], ["sigma_hk", "sigma_incl", "sigma_incl_near_mmr"]]
active_params_box = [[(log(1.), log(5.)), (log(1.), log(5.)), (3., 20.), (-0.5, 1.5), (0., 0.3)], [(0.5*ExoplanetsSysSim.earth_radius, 10.*ExoplanetsSysSim.earth_radius),(1., 4.), (-6., 0.), (-6., 0.), (0.1, 1.0)], [(0., 0.1), (0., 5.), (0., 5.)]] #search ranges for all of the active parameters
    #[[(log(1.), log(5.)), (log(1.), log(5.)), (3., 20.), (-0.5, 1.5), (0., 0.3)], [(0.5*ExoplanetsSysSim.earth_radius, 10.*ExoplanetsSysSim.earth_radius),(1., 4.), (-6., 0.), (-6., 0.), (0.1, 1.0)], [(0., 0.1), (0., 5.), (0., 5.)]] #search ranges for all of the active parameters

active_param_draws = [Float64[], Float64[], Float64[]] #[Float64[], Float64[], Float64[]]
for (i,active_param_block) in enumerate(active_param_keys)
    for (j,param_key) in enumerate(active_param_block)
        active_param_draw = active_params_box[i][j][1] + (active_params_box[i][j][2] - active_params_box[i][j][1])*rand(1)
        append!(active_param_draws[i], active_param_draw)
        add_param_fixed(sim_param, param_key, active_param_draw[1])
    end
end

println("# All active parameters: ", active_param_keys)
println(f, "# All active parameters: ", active_param_keys)
println(f, "# Starting active parameter values: ", active_param_draws)
println(f, "# Optimization active parameters search bounds: ", active_params_box)
println(f, "# Method: adaptive_de_rand_1_bin_radiuslimited")
println(f, "# PopulationSize: ", Pop_per_param, " (per param)")
println(f, "# Format: Active_params: [active parameter values]")
println(f, "# Format: Dist: [distances][total distance]")
println(f, "# Format: Dist_weighted: [weighted distances][total weighted distance]")
println(f, "# Distances used: ", use_KS_or_AD)
println(f, "# AD_mod: ", AD_mod)

#To start the block-optimization:

# Pkg.add("BlackBoxOptim")       # only need to do these once
# Pkg.checkout("BlackBoxOptim")  # needed to get the lastest version
using BlackBoxOptim              # see https://github.com/robertfeldt/BlackBoxOptim.jl for documentation

tic()
for c in 1:cycles
    println("# Cycle: ", c)
    println(f, "#")
    println(f, "# Cycle: ", c)
    println(f, "#")

    tic()
    for (i,active_param_block) in enumerate(active_param_keys)
        println("# Active parameters: ", active_param_block)
        println(f, "# Active parameters: ", active_param_block)

        #To make the current block of parameters active:
        for (j,param_key) in enumerate(active_param_block)
            add_param_active(sim_param, param_key, active_param_draws[i][j])
        end

        tic()
        opt_result = bboptimize(active_params -> target_function(active_params, use_KS_or_AD, Kep_or_Sim ; AD_mod=AD_mod, weights=weights, all_dist=false, save_dist=true); SearchRange = active_params_box[i], NumDimensions = length(active_param_block), Method = :adaptive_de_rand_1_bin_radiuslimited, PopulationSize = length(active_param_block)*Pop_per_param, MaxFuncEvals = max_evals, TargetFitness = target_fitness, FitnessTolerance = target_fitness_std, TraceMode = :verbose)
        t_block = toc()
        active_params_best = best_candidate(opt_result)
        best_distance = best_fitness(opt_result)

        println(f, "# Best parameters block: ", active_params_best)
        println(f, "# Best fitness block: ", best_distance)
        println(f, "# Block elapsed time: ", t_block, " seconds")
        println(f, "#")

        #To make the block of parameters fixed at the best values found:
        for (j,param_key) in enumerate(active_param_block)
            add_param_fixed(sim_param, param_key, active_params_best[j])
        end
    end
    t_cycle = toc()

    println(f, "# Cycle elapsed time: ", t_cycle, " seconds")
end
t_optim = toc()

println(f, "#")
println(f, "# Total elapsed time: ", t_optim, " seconds")
close(f)