text:
include("../src/hypergraph-affinity-functions.jl")

dataset = "congress-bills"
r = 25
H, classes = read_hypergraph_data(dataset,r)
# 1 = democrat
# 2 = republican
m,n = size(H)
classes = classes .- 1
# 0 = democrat = 2 mod 0
# 1 = republican

alpha, B1, B2, H1, H2, R1, R2, N = Hypergraph_to_Scores(H,classes,25)
class1 = "Republican"
class2 = "Democrat"
hr, hd = clique_expansion_homophily(N,r)

## Plot
lw = 2.5
ms = 9
gfs = 15
tfs = 15
titlesize = 14
s1 = 350
s2 = 300
k = 5
color1 = :red
color2 = :blue
title = "$k-person bills"
title = ""
leg = false
xlab = "Affinity type t"
ylab = "Affinity"
p = Plots.plot(fontfamily = "helvetica",linewidth = 2,legend = false,grid = false, size = (s1,s2),color = :gray)
title!(p,title,titlefontsize = titlesize)
xaxis!(p,xlab, xticks = k/5:k/5:k,tickfontsize=tfs,guidefontsize = gfs)
yaxis!(p,ylab,tickfontsize=tfs,guidefontsize = gfs)

plot!(p,1:k,H2[k,1:k],color = color2, linewidth = lw, markersize = ms,
    markershape = :circle, markerstrokecolor = color2, label = class2)
plot!(p,1:k,H1[k,1:k],color = color1, linewidth = lw, markersize = ms,
    markershape = :circle,markerstrokecolor = color1, label = class1)

savefig("Figures/Affinities_$(k)_Congress.pdf")

## Next
k = 15
s1 = 400
s2 = 300
title = "$k-person bills"
title = ""
leg = false
xlab = "Affinity type t"
ylab = "Affinity / Baseline"
p = Plots.plot([1; k],[1; 1],yscale = :log10,fontfamily = "helvetica",linewidth = 2,legend = false,grid = false, size = (s1,s2),color = :gray)
title!(p,title,titlefontsize = titlesize)
xaxis!(p,xlab, xticks = k/5:k/5:k,tickfontsize=tfs,guidefontsize = gfs)
yaxis!(p,ylab,tickfontsize=tfs,guidefontsize = gfs)
plot!(p,1:k,R1[k,1:k],color = color1, linewidth = lw, markersize = ms,
    markershape = :circle,markerstrokecolor = color1, label = class1)
plot!(p,1:k,R2[k,1:k],color = color2, linewidth = lw, markersize = ms,
    markershape = :circle, markerstrokecolor = color2, label = class2)

savefig("Figures/Ratios_$(k)_Congress.pdf")

## Group Homophily Index, print usable LaTeX code for building table
S1,S2  = GHI(H,classes,r)
print("\n \$k\$ ")
start = 5
for k = start:20
    print("& $(k) ")
end
print("\\\\\n \\midrule \nRep. GHI ")
for k = start:20
    print("& $(round(Int64,S1[k])) ")
end
print("\\\\ \n Dem. GHI")
for k = start:20
    print("& $(round(Int64,S2[k])) ")
end
print("\\\\ \n")

max_stars_repo_name: nveldt/HypergraphHomophily

max_issues_repo_name: nveldt/HypergraphHomophily

max_forks_repo_name: nveldt/HypergraphHomophily


text:
"""
# Example:
     find_in = "D:\\GoogleDrive\\1Programming\\DATA"
     regexppattern = r".*\\.mat"
     find = "both"
     method = "readdir"
     O = datalist(regexppattern, find_in, method, find);
# Output:
     O.fullpath is a string vector of the paths of the desired files or folders.
"""
function datalist(regexppattern::Regex, find_in::String,
                  method::String="readdir", find::String="both")
# CHECKPOINT:
# for the setting:
#     find_in = "D:\\GoogleDrive\\1Programming\\DATA"
#     find = "both"
#     method = "walkdir"
# The following results are verified as identical in their number:
#     datalist(r".*\.mat", find_in, method, find); # julia
#     datalist('*.mat',find_in,'Search','**'); % in matlab
#
#     datalist(r"\[.*\].*\.mat", find_in, method, find);
#     datalist('*[*]*.mat',find_in,'Search','**')

# code start here
    findfolders = true;
    findfiles = true;
    if find == "folder"
        findfiles = false;
    elseif find == "file"
        findfolders = false;
    end

    if method == "readdir";
        result_paths = readdir(find_in;join=true);
    elseif method == "walkdir"
        result_iterator = walkdir(find_in);
        result_paths = Vector{String}(); # or String[]
        for (root, dirs, files) in result_iterator
            # println("Directories in $root")
            in_folder_i = Vector{String}();
            if findfolders
                for dir in dirs
                    # path to directories
                    in_folder_i = vcat(in_folder_i,joinpath(root, dir));
                    # result_paths = push!(result_paths, joinpath(root, dir));
                end
            end
            if findfiles
                # println("Files in $root")
                for file in files
                    in_folder_i = vcat(in_folder_i,joinpath(root, file));
                    # println(joinpath(root, file)) # path to files
                end

            end
            result_paths = vcat(result_paths,in_folder_i);
        end
    end
    # check if the file/folder names match the regular expression pattern, and
    # return a 1-d true false array.
    desired_ind = occursin.(regexppattern,basename.(result_paths));
    O = (fullpath = sort(result_paths[desired_ind]),);
    return O
end

# test:
# find_in = "D:\\GoogleDrive\\1Programming\\DATA"
# regexppattern = r"\[.*\].*\.mat"
# find = "both"
# method = "walkdir"
# O = datalist(regexppattern, find_in, method, find).fullpath;
# lenpaths = length(O)
#
# for i = 1:lenpaths
#     println(O[i]);
# end

ext: jl


text:
using EpdTest
using Test

@testset "EpdTest.jl" begin
    μ,σ,p = 0., 1., 1.
    y = [0.2, -1., 2.2, 3.1]
    @test EpdTest.epdTest(y, μ, σ, p) === -1.4366494480275978
end


text:
print("Please enter sleep duration in seconds: ")
input = int(readline(STDIN))
println("Sleeping...")
sleep(input)
println("Awake!")

max_stars_repo_name: ethansaxenian/RosettaDecode

max_issues_repo_name: ethansaxenian/RosettaDecode

max_forks_repo_name: ethansaxenian/RosettaDecode


text:
# A box advection test to visualise how different filters work

using MPI
using OrderedCollections
using Plots
using StaticArrays
using Printf

using CLIMAParameters
struct EarthParameterSet <: AbstractEarthParameterSet end
const param_set = EarthParameterSet()

using ClimateMachine
using ClimateMachine.Mesh.Topologies
using ClimateMachine.Mesh.Grids
using ClimateMachine.DGMethods
using ClimateMachine.DGMethods.NumericalFluxes
using ClimateMachine.BalanceLaws:
    BalanceLaw, Prognostic, Auxiliary, Gradient, GradientFlux
using ClimateMachine.Mesh.Geometry: LocalGeometry
using ClimateMachine.Mesh.Filters
using ClimateMachine.MPIStateArrays
using ClimateMachine.GenericCallbacks
using ClimateMachine.ODESolvers
using ClimateMachine.VariableTemplates
using ClimateMachine.SingleStackUtils

import ClimateMachine.BalanceLaws:
    vars_state,
    source!,
    flux_second_order!,
    flux_first_order!,
    compute_gradient_argument!,
    compute_gradient_flux!,
    update_auxiliary_state!,
    nodal_init_state_auxiliary!,
    init_state_prognostic!,
    boundary_conditions,
    wavespeed

ClimateMachine.init(; disable_gpu = true, log_level = "warn");
const clima_dir = dirname(dirname(pathof(ClimateMachine)));
include(joinpath(clima_dir, "docs", "plothelpers.jl"));

Base.@kwdef struct Box1D{FT, _init_q, _amplitude, _velo} <: BalanceLaw
    param_set::AbstractParameterSet = param_set
    init_q::FT = _init_q
    amplitude::FT = _amplitude
    velo::FT = _velo
end

vars_state(::Box1D, ::Auxiliary, FT) = @vars(z_dim::FT);
vars_state(::Box1D, ::Prognostic, FT) = @vars(q::FT);
vars_state(::Box1D, ::Gradient, FT) = @vars();
vars_state(::Box1D, ::GradientFlux, FT) = @vars();

function wavespeed(
    ::Box1D{FT, _init_q, _amplitude, _velo},
    _...,
) where {FT, _init_q, _amplitude, _velo}
    return _velo
end

function nodal_init_state_auxiliary!(
    m::Box1D,
    aux::Vars,
    tmp::Vars,
    geom::LocalGeometry,
)
    aux.z_dim = geom.coord[3]
end;

function init_state_prognostic!(
    m::Box1D,
    state::Vars,
    aux::Vars,
    localgeo,
    t::Real,
)
    if aux.z_dim >= 75 && aux.z_dim <= 125
        state.q = m.init_q + m.amplitude
    else
        state.q = m.init_q
    end
end;

function update_auxiliary_state!(
    dg::DGModel,
    m::Box1D,
    Q::MPIStateArray,
    t::Real,
    elems::UnitRange,
)
    return true
end;

function source!(m::Box1D, _...) end;

@inline function flux_first_order!(
    m::Box1D,
    flux::Grad,
    state::Vars,
    aux::Vars,
    t::Real,
    _...,
)
    FT = eltype(state)
    @inbounds begin
        flux.q = SVector(FT(0), FT(0), state.q * m.velo)
    end
end

@inline function flux_second_order!(
    m::Box1D,
    flux::Grad,
    state::Vars,
    diffusive::Vars,
    hyperdiffusive::Vars,
    aux::Vars,
    t::Real,
) end

@inline function flux_second_order!(
    m::Box1D,
    flux::Grad,
    state::Vars,
    τ,
    d_h_tot,
) end

boundary_condtions(m::Box1D) = ()

function run_box1D(
    N_poly::Int,
    init_q::FT,
    amplitude::FT,
    velo::FT,
    plot_name::String;
    tmar_filter::Bool = false,
    cutoff_filter::Bool = false,
    exp_filter::Bool = false,
    boyd_filter::Bool = false,
    cutoff_param::Int = 1,
    exp_param_1::Int = 0,
    exp_param_2::Int = 32,
    boyd_param_1::Int = 0,
    boyd_param_2::Int = 32,
    numerical_flux_first_order = CentralNumericalFluxFirstOrder(),
) where {FT}
    N_poly = N_poly
    nelem = 128
    zmax = FT(350)

    m = Box1D{FT, init_q, amplitude, velo}()

    driver_config = ClimateMachine.SingleStackConfiguration(
        "Box1D",
        N_poly,
        nelem,
        zmax,
        param_set,
        m,
        numerical_flux_first_order = numerical_flux_first_order,
        boundary = ((0, 0), (0, 0), (0, 0)),
        periodicity = (true, true, true),
    )

    t0 = FT(0)
    timeend = FT(450)

    Δ = min_node_distance(driver_config.grid, VerticalDirection())
    max_vel = m.velo
    dt = Δ / max_vel

    solver_config = ClimateMachine.SolverConfiguration(
        t0,
        timeend,
        driver_config,
        ode_dt = dt,
    )
    grid = solver_config.dg.grid
    Q = solver_config.Q
    aux = solver_config.dg.state_auxiliary

    output_dir = @__DIR__
    mkpath(output_dir)

    z_label = "z"
    z = get_z(grid)

    # store initial condition at ``t=0``
    dons_arr = Dict[dict_of_nodal_states(solver_config)]
    time_data = FT[0]                                      # store time data

    # output
    output_freq = floor(Int, timeend / dt) + 10

    cb_output = GenericCallbacks.EveryXSimulationSteps(output_freq) do
        push!(dons_arr, dict_of_nodal_states(solver_config))
        push!(time_data, gettime(solver_config.solver))
        nothing
    end

    filter_freq = 1
    # tmar filter
    cb_tmar =
        GenericCallbacks.EveryXSimulationSteps(filter_freq) do (init = false)
            Filters.apply!(
                solver_config.Q,
                (:q,),
                solver_config.dg.grid,
                TMARFilter(),
            )
            nothing
        end
    # cutoff filter
    cb_cutoff =
        GenericCallbacks.EveryXSimulationSteps(filter_freq) do (init = false)
            Filters.apply!(
                solver_config.Q,
                (:q,),
                solver_config.dg.grid,
                CutoffFilter(solver_config.dg.grid, cutoff_param),
            )
            nothing
        end
    # exponential filter
    cb_exp =
        GenericCallbacks.EveryXSimulationSteps(filter_freq) do (init = false)
            Filters.apply!(
                solver_config.Q,
                (:q,),
                solver_config.dg.grid,
                ExponentialFilter(
                    solver_config.dg.grid,
                    exp_param_1,
                    exp_param_2,
                ),
            )
            nothing
        end
    # exponential filter
    cb_boyd =
        GenericCallbacks.EveryXSimulationSteps(filter_freq) do (init = false)
            Filters.apply!(
                solver_config.Q,
                (:q,),
                solver_config.dg.grid,
                BoydVandevenFilter(
                    solver_config.dg.grid,
                    boyd_param_1,
                    boyd_param_2,
                ),
            )
            nothing
        end

    user_cb_arr = [cb_output]
    if tmar_filter
        push!(user_cb_arr, cb_tmar)
    end
    if cutoff_filter
        push!(user_cb_arr, cb_cutoff)
    end
    if exp_filter
        push!(user_cb_arr, cb_exp)
    end
    if boyd_filter
        push!(user_cb_arr, cb_boyd)
    end
    user_cb = (user_cb_arr...,)

    initial_mass = weightedsum(solver_config.Q)
    ClimateMachine.invoke!(solver_config; user_callbacks = (user_cb))
    final_mass = weightedsum(solver_config.Q)
    @info @sprintf(
        """
Mass Conservation:
    initial mass          = %.16e
    final mass            = %.16e
    difference            = %.16e
    normalized difference = %.16e""",
        initial_mass,
        final_mass,
        final_mass - initial_mass,
        (final_mass - initial_mass) / initial_mass
    )

    push!(dons_arr, dict_of_nodal_states(solver_config))
    push!(time_data, gettime(solver_config.solver))

    export_plot(
        z,
        time_data,
        dons_arr,
        ("q",),
        joinpath(output_dir, plot_name);
        xlabel = "x",
        ylabel = "q",
        horiz_layout = true,
    )

end

max_forks_repo_name: jm-c/CLIMA


text:
function getRowKernel{T <: FloatingPoint}(::Type{T})
    nType = T == Float64 ? "double" : "float"

    return "
        #if defined(cl_khr_fp64)  // Khronos extension available?
        #pragma OPENCL EXTENSION cl_khr_fp64 : enable
        #elif defined(cl_amd_fp64)  // AMD extension available?
        #pragma OPENCL EXTENSION cl_amd_fp64 : enable
        #endif
        #define number $nType

        #define M(x,y) a[y*D2 + x]
        #define A(x,y) b[y*D2 + x]
        #define F(x,y) c[y*D2 + x]
        #define W(x,y) d[y*D2 + x]

        #define Out(x,y) out[y*D2 + x]

        __kernel void calcRow(
                      __global const number *a,
                      __global const number *b,
                      __global const number *c,
                      __global const number *d,
                      __global number *out,
                      const int D1,
                      const int D2,
                      const number m,
                      const number _a,
                      const number f,
                      const number w) {

        int i = get_global_id(0);
        int j = get_global_id(1);

        Out(i, j) = pow(M(i, j), m) * pow(A(i, j), _a) * pow(F(i, j), f) * pow(W(i, j), w);
    }
"
end

function calcRowCL!{T <: FloatingPoint}(
    a_buff, b_buff, c_buff, d_buff,
    out_buff,
    a :: Real, b :: Real, c :: Real, d :: Real, d1 :: Int64, d2 :: Int64,
    ctx, queue, program, ::Type{T})

    k = Kernel(program, "calcRow")

    call(queue, k, (d1,d2), nothing, a_buff, b_buff, c_buff, d_buff, out_buff, int32(d1), int32(d2), convert(T, a), convert(T, b), convert(T, c), convert(T, d))
end


text:
#
# types.jl --
#
# Type definitions for XPA package.
#
#------------------------------------------------------------------------------
#
# This file is part of XPA.jl released under the MIT "expat" license.
# Copyright (C) 2016-2019, Éric Thiébaut (https://github.com/emmt/XPA.jl).
#

const Byte = UInt8
const NULL = Ptr{Byte}(0)

"""

`XPA.SUCCESS` and `XPA.FAILURE` are the possible values returned by the
callbacks of an XPA server.

"""
const SUCCESS = convert(Cint,  0)
const FAILURE = convert(Cint, -1)
@doc @doc(SUCCESS) FAILURE

# Server mode flags for receive, send, info.  These constants are defined
# in the *private* header file `xpap.h`.
const MODE_BUF     = 1 # XPA_MODE_BUF in `xpap.h`
const MODE_FILLBUF = 2 # XPA_MODE_FILLBUF in `xpap.h`
const MODE_FREEBUF = 4 # XPA_MODE_FREEBUF in `xpap.h`
const MODE_ACL     = 8 # XPA_MODE_ACL in `xpap.h`

"""

Abstract type `XPA.Handle` is the super type of client ([`XPA.Client`](@ref))
and server ([`XPA.Server`](@ref)) connections in the XPA Messaging System.

"""
abstract type Handle end

# All concrete types derived from `Handle` have a `ptr` field which is NULL if
# handle has been closed (or not yet open).
Base.isopen(xpa::Handle) = xpa.ptr != C_NULL

"""

An instance of the mutable structure `XPA.Client` represents a client
connection in the XPA Messaging System.

"""
mutable struct Client <: Handle # must be mutable to be finalized
    ptr::Ptr{Cvoid} # pointer to XPARec structure
end

"""

`XPA.TupleOf{T}` represents a tuple of any number of elements of type `T`, it
is an alias for `Tuple{Vararg{T}}`

"""
const TupleOf{T} = Tuple{Vararg{T}}

"""

`XPA.Reply` is used to store the answer(s) of [`XPA.get`](@ref) and
[`XPA.set`](@ref) requests.  Method `length` applied to an object of type
`Reply` yields the number of replies.  Methods [`XPA.get_data`](@ref),
[`XPA.get_server`](@ref) and [`XPA.get_message`](@ref) can be used to retrieve
the contents of an object of type `XPA.Reply`.

"""
mutable struct Reply
    replies::Int
    lengths::Vector{Csize_t}
    buffers::Vector{Ptr{Byte}}
end

Base.length(rep::Reply) = rep.replies

"""

An instance of the mutable structure `XPA.Server` represents a server
connection in the XPA Messaging System.

"""
mutable struct Server <: Handle # must be mutable to be finalized
    ptr::Ptr{Cvoid} # pointer to XPARec structure
end

abstract type Callback end

"""

An instance of the `XPA.SendCallback` structure represents a callback called to
serve an [`XPA.get`](@ref) request.

"""
mutable struct SendCallback{T,F<:Function} <: Callback
    # must be mutable because pointer_from_objref is used to recover it
    send::F        # function to call on `XPA.get` requests
    data::T        # client data
    acl::Bool      # enable access control
end

"""

An instance of the `XPA.SendBuffer` structure is provided to send callbacks to
record the addresses where to store the address and size of the data associated
to the answer of an [`XPA.get`](@ref) request.  A send callback shall use
[`XPA.store!`](@ref) to set the buffer contents.

See also [`XPA.store!`](@ref), [`XPA.get`](@ref), [`XPA.Server`](@ref)
and [`XPA.SendCallback`](@ref).

"""
struct SendBuffer
    bufptr::Ptr{Ptr{Byte}}
    lenptr::Ptr{Csize_t}
end

"""

An instance of the `XPA.ReceiveCallback` structure represents a callback called
to serve an [`XPA.set`](@ref) request.

"""
mutable struct ReceiveCallback{T,F<:Function} <: Callback
    # must be mutable because pointer_from_objref is used to recover it
    recv::F        # function to call on `XPA.set` requests
    data::T        # client data
    acl::Bool      # enable access control
end

"""

An instance of the `XPA.ReceiveBuffer` structure is provided to receive
callbacks to record the address and the size of the data sent by an
[`XPA.set`](@ref) request.  Methods `pointer(buf)` and `sizeof(buf)` can be
used to query the address and the number of bytes of the buffer `buf`.

See also [`XPA.get`](@ref), [`XPA.Server`](@ref) and
[`XPA.ReceiveCallback`](@ref).

"""
struct ReceiveBuffer
    ptr::Ptr{Byte}
    len::Int
    # Inner constructor that checks parameters at construction time.
    function ReceiveBuffer(ptr::Ptr{Byte}, len::Integer)
        @assert ptr == NULL ? len == 0 : len ≥ 0
        return new(ptr, len)
    end
end
Base.sizeof(buf::ReceiveBuffer)::Int = buf.len
Base.pointer(buf::ReceiveBuffer) = buf.ptr

"""

An instance of the `XPA.AccessPoint` structure represents an available XPA
server.  A vector of such instances is returned by the [`XPA.list`](@ref)
utility.

"""
struct AccessPoint
    class::String # class of the access point
    name::String  # name of the access point
    addr::String  # socket access method (host:port for inet,
                  # file for local/unix)
    user::String  # user name of access point owner
    access::UInt  # allowed access
end

"""

`XPA.NullBuffer` is a singleton type representing a NULL-buffer when sending
data to a server.

"""
struct NullBuffer end


text:
module DocumenterGitHubActionsExample

# Write your package code here.
function foo(x)
    y = x^2
    if x % 2 == 0
        y = -1
    end
    return y+1
end

end

max_stars_repo_name: mortenpi/DocumenterGitHubActionsExample

max_issues_repo_name: mortenpi/DocumenterGitHubActionsExample

max_forks_repo_name: mortenpi/DocumenterGitHubActionsExample


text:
const registry=Dict{Tuple, Any}()
const refs=Set()  # Collection of darray identities created on this node

let DID::Int = 1
    global next_did
    next_did() = (id = DID; DID += 1; (myid(), id))
end

"""
    next_did()

Produces an incrementing ID that will be used for DArrays.
"""
next_did

release_localpart(id::Tuple) = (delete!(registry, id); nothing)
release_localpart(d) = release_localpart(d.id)

function close_by_id(id, pids)
#   @async println("Finalizer for : ", id)
    global refs
    @sync begin
        for p in pids
            @async remotecall_fetch(release_localpart, p, id)
        end
        if !(myid() in pids)
            release_localpart(id)
        end
    end
    delete!(refs, id)
    nothing
end

function Base.close(d::DArray)
#    @async println("close : ", d.id, ", object_id : ", object_id(d), ", myid : ", myid() )
    if (myid() == d.id[1]) && d.release
        @async close_by_id(d.id, d.pids)
        d.release = false
    end
    nothing
end

function d_closeall()
    crefs = copy(refs)
    for id in crefs
        if id[1] ==  myid() # sanity check
            haskey(registry, id) && close(d_from_weakref_or_d(id))
            yield()
        end
    end
end

"""
    procs(d::DArray)

Get the vector of processes storing pieces of DArray `d`.
"""
Distributed.procs(d::DArray)    = d.pids
Distributed.procs(d::SubDArray) = procs(parent(d))

"""
    localpart(A)

The identity when input is not distributed
"""
localpart(A) = A


text:
using LinearAlgebra
@testset "rational" begin
    A = [-3 1; 5 6.6]

    # Simple calls with coefficients
    den_coeff = [3; -1; 2.0; 0.1; 1.2; 0.9; -0.3; 1.6]
    num_coeff = [1.4, 2.1, -4, -0.1, 0.2, -1.2, 2.1, 0.4, 1.4]
    den =
        den_coeff[1] * I +
        den_coeff[2] * A +
        den_coeff[3] * A^2 +
        den_coeff[4] * A^3 +
        den_coeff[5] * A^4 +
        den_coeff[6] * A^5 +
        den_coeff[7] * A^6 +
        den_coeff[8] * A^7
    num =
        num_coeff[1] * I +
        num_coeff[2] * A +
        num_coeff[3] * A^2 +
        num_coeff[4] * A^3 +
        num_coeff[5] * A^4 +
        num_coeff[6] * A^5 +
        num_coeff[7] * A^6 +
        num_coeff[8] * A^7 +
        num_coeff[9] * A^8
    P = num \ den

    (graph, cref) = graph_rational(den_coeff, num_coeff, graph_monomial)
    @test eval_graph(graph, A) ≈ P

    (graph, cref) = graph_rational(den_coeff, num_coeff, graph_ps)
    @test eval_graph(graph, A) ≈ P

    # Generate more complicated graphs before
    x1a = [1.2; 4.4]
    x1b = [0.3; -0.3]
    x2a = [2.1; 0.1; -0.1]
    x2b = [0.3; 0.2; 0.1]
    x = [(x1a, x1b), (x2a, x2b)]
    z = [1.0; 0.1; -0.1; 0.01]
    den_coeffs = [(x1a, x1b), (x2a, x2b)]
    β1 = (x1a[1] * I + x1a[2] * A) * (x1b[1] * I + x1b[2] * A)
    β2 =
        (x2a[1] * I + x2a[2] * A + x2a[3] * β1) *
        (x2b[1] * I + x2b[2] * A + x2b[3] * β1)
    den = z[1] * I + z[2] * A + z[3] * β1 + z[4] * β2
    (den_graph, den_cref) = graph_degopt(x, z)

    x1a = [2.1; 1.2]
    x1b = [1.3; 0.3]
    x2a = [0.1; -1; 1]
    x2b = [0; 2; -0.3]
    x = [(x1a, x1b), (x2a, x2b)]
    z = [1.0; -1; 0.1; -0.01]
    num_coeffs = [(x1a, x1b), (x2a, x2b)]
    β1 = (x1a[1] * I + x1a[2] * A) * (x1b[1] * I + x1b[2] * A)
    β2 =
        (x2a[1] * I + x2a[2] * A + x2a[3] * β1) *
        (x2b[1] * I + x2b[2] * A + x2b[3] * β1)
    num = z[1] * I + z[2] * A + z[3] * β1 + z[4] * β2
    (num_graph, num_cref) = graph_degopt(x, z)

    P = num \ den
    (graph, cref) = graph_rational(
        den_graph,
        num_graph,
        den_cref = den_cref,
        num_cref = num_cref,
    )
    @test eval_graph(graph, A) ≈ P
end


text:
const TRACE = haskey(ENV, "TRACE")
"Display a trace message. Only results in actual printing if the TRACE environment variable
is set."
@inline function trace(io::IO, msg...; prefix="TRACE: ", line=true)
    @static if TRACE
        Base.print_with_color(:cyan, io, prefix, chomp(string(msg...)))
        if line
            println(io)
        end
    end
end
@inline trace(msg...; kwargs...) = trace(STDERR, msg...; kwargs...)

const DEBUG = TRACE || haskey(ENV, "DEBUG")
"Display a debug message. Only results in actual printing if the TRACE or DEBUG environment
variable is set."
@inline function debug(io::IO, msg...; prefix="DEBUG: ", line=true)
    @static if DEBUG
        Base.print_with_color(:green, io, prefix, chomp(string(msg...)))
        if line
            println(io)
        end
    end
end
@inline debug(msg...; kwargs...) = debug(STDERR, msg...; kwargs...)

"Create an indented string from any value (instead of escaping endlines as \n)"
function repr_indented(ex; prefix=" "^7, abbrev=true)
    io = IOBuffer()
    print(io, ex)
    str = takebuf_string(io)

    # Limit output
    if abbrev && length(str) > 256
        if isa(ex, Array)
            T = eltype(ex)
            dims = join(size(ex), " by ")
            if isbits(T) && zeros(ex) == ex     # NOTE: not all types have zero() defined
                str = "$T[$dims zeros]"
            else
                str = "$T[$dims elements]"
            end
        else
            if contains(strip(str), "\n")
                str = str[1:100] * "…\n\n[snip]\n\n…" * str[end-100:end]
            else
                str = str[1:100] * "…" * str[end-100:end]
            end
        end
    end

    lines = split(strip(str), '\n')
    if length(lines) > 1
        for i = 1:length(lines)
            lines[i] = prefix * lines[i]
        end

        lines[1] = "\"\n" * lines[1]
        lines[length(lines)] = lines[length(lines)] * "\""

        return join(lines, '\n')
    else
        return str
    end
end

function __init_logging__()
    # TODO: assign TRACE and DEBUG at run-time, not using the pre-compiled code
    #       when the values are different?
    #       or make it work like CPU_CORES dose after Julia/#16219

    if TRACE
        trace("CUDAdrv.jl is running in trace mode, this will generate a lot of additional output")
    elseif DEBUG
        debug("CUDAdrv.jl is running in debug mode, this will generate additional output")
        debug("Run with TRACE=1 to enable even more output")
    end
end


text:
using Distributions, UncertainData, Combinatorics

# Test all combinations of different types of uncertain values
M = MixtureModel([Normal(3, 0.2), Normal(2, 1)])

r1 = UncertainValue(Normal, rand(), rand())
r2 = UncertainValue(rand(M, 10000))
r3 = UncertainValue(Normal, rand(Normal(4, 3.2), 10000))
r4 = CertainValue(2.2)
r5 = CertainValue(2)

uvals = [r1; r2; r3]

# Find all the different ways of combining two uncertain values in `uvals`.
# We let the order matter, so that if the order of operations for some reason causes
# things to fail when using different types of uncertain values, the errors are caught.
perms = permutations(1:3, 2) |> collect

# A random number
x = rand()

# Number of draws when deviating from default n = 10000
n = 10

for perm in perms
    r1, r2 = uvals[perm[1]], uvals[perm[2]]

    # Addition
    @test r1 + r2 isa AbstractUncertainValue
    @test x + r2 isa AbstractUncertainValue
    @test r1 + x isa AbstractUncertainValue
    @test +(r1, r2, n) isa AbstractUncertainValue
    @test +(x, r2, n) isa AbstractUncertainValue
    @test +(r1, x, n) isa AbstractUncertainValue

    # Subtraction
    @test r1 - r2 isa AbstractUncertainValue
    @test x - r2 isa AbstractUncertainValue
    @test r1 - x isa AbstractUncertainValue
    @test -(r1, r2, n) isa AbstractUncertainValue
    @test -(x, r2, n) isa AbstractUncertainValue
    @test -(r1, x, n) isa AbstractUncertainValue

    # Multiplication
    @test r1 * r2 isa AbstractUncertainValue
    @test x * r2 isa AbstractUncertainValue
    @test r1 * x isa AbstractUncertainValue
    @test *(r1, r2, n) isa AbstractUncertainValue
    @test *(x, r2, n) isa AbstractUncertainValue
    @test *(r1, x, n) isa AbstractUncertainValue

    # Division
    @test r1 / r2 isa AbstractUncertainValue
    @test x / r2 isa AbstractUncertainValue
    @test r1 / x isa AbstractUncertainValue
    @test /(r1, r2, n) isa AbstractUncertainValue
    @test /(x, r2, n) isa AbstractUncertainValue
    @test /(r1, x, n) isa AbstractUncertainValue

    # Exponentiation
    #@test r1 ^ r2 isa AbstractUncertainValue # implement for complex numbers.
end

# Elementary operations when one or both arguments is a vector of uncertain values
for val in uvals
    @test val + uvals isa Vector{<:AbstractUncertainValue}
    @test uvals + val isa Vector{<:AbstractUncertainValue}
    @test uvals + uvals isa Vector{<:AbstractUncertainValue}
end


text:
function RCI(name = :RCI)
    model = BondGraph(name)
    C = Component(:C)
    R = Component(:R)
    I = Component(:I)
    SS = Component(:SS)
    zero_law = EqualEffort()

    add_node!(model, [C, R, I, SS, zero_law])
    connect!(model, R, zero_law)
    connect!(model, C, zero_law)
    connect!(model, zero_law, I)
    connect!(model, zero_law, SS)

    model
end

# Based on https://bondgraphtools.readthedocs.io/en/latest/tutorials/RC.html
@testset "BondGraph Construction" begin
    model = BondGraph(:RC)
    C = Component(:C)
    R = Component(:R)
    zero_law = EqualEffort()

    add_node!(model, [R, C, zero_law])
    @test R in model.nodes
    @test C in model.nodes
    @test zero_law in model.nodes

    b1 = connect!(model, R, zero_law)
    b2 = connect!(model, zero_law, C)
    @test b1 in model.bonds
    @test b2 in model.bonds
end

@testset "BondGraph Modification" begin
    model = BondGraph(:RCI)
    C = Component(:C)
    R = Component(:R)
    I = Component(:I)
    SS = Component(:SS)
    zero_law = EqualEffort()
    one_law = EqualFlow()

    add_node!(model, [C, R, I, SS, zero_law, one_law])
    remove_node!(model, [SS, one_law])
    @test !(SS in model.nodes)
    @test !(one_law in model.nodes)

    connect!(model, R, zero_law)
    connect!(model, C, zero_law)

    @test I.freeports == [true]
    b1 = connect!(model, zero_law, I)
    @test ne(model) == 3
    @test b1 in model.bonds
    @test I.freeports == [false]

    disconnect!(model, I, zero_law) # test disconnect when node order is swapped
    @test ne(model) == 2
    @test !(b1 in model.bonds)
    @test I.freeports == [true]

    connect!(model, zero_law, I)
    swap!(model, zero_law, one_law)
    @test I.freeports == [false]
    @test one_law in model.nodes
    @test inneighbors(model, one_law) == [R, C]
    @test outneighbors(model, one_law) == [I]
end

@testset "Construction Failure" begin
    model = BondGraph(:RC)
    C = Component(:C)
    R = Component(:R)
    zero_law = EqualEffort()

    add_node!(model, [R, C, zero_law])
    @test_logs (:warn, "R:R already in model") add_node!(model, R)
    @test_logs (:warn, "0 already in model") add_node!(model, zero_law)

    bond = connect!(model, R, zero_law)
    @test_throws ErrorException connect!(model, R, zero_law)
    @test_throws ErrorException connect!(model, C, R)

    one_law = EqualFlow()
    @test_logs (:warn, "1 not in model") remove_node!(model, one_law)

    tf = Component(:TF)
    add_node!(model, tf)
    @test_throws ErrorException swap!(model, tf, C)

    # if inserting a node fails, the original nodes should still remain connected
    @test_throws ErrorException insert_node!(model, bond, Component(:I))
    @test has_edge(model, bond)
end

@testset "Chemical reaction" begin
    model = BondGraph(:Chemical)
    A = Component(:C, :A)
    B = Component(:C, :B)
    C = Component(:C, :C)
    D = Component(:C, :D)
    Re = Component(:Re, :Reaction, numports = 2)
    J_AB = EqualFlow()
    J_CD = EqualFlow()

    add_node!(model, [A, B, C, D, Re, J_AB, J_CD])
    connect!(model, A, J_AB)
    connect!(model, B, J_AB)
    connect!(model, C, J_CD)
    connect!(model, D, J_CD)

    @test freeports(Re) == [true, true]
    @test freeports(J_AB) == [false, false]

    # Connecting junctions to specific ports in Re
    connect!(model, Re, J_CD, srcportindex = 2)
    @test freeports(Re) == [true, false]

    # connecting to a full port should fail
    @test_throws ErrorException connect!(model, J_AB, Re, dstportindex = 2)

    connect!(model, J_AB, Re, dstportindex = 1)
    @test freeports(Re) == [false, false]

    @test nv(model) == 7
    @test ne(model) == 6
end

@testset "Standard components" begin
    tf = Component(:TF, :n)
    @test tf isa Component{2}
    @test tf.type == :TF
    @test numports(tf) == 2

    r = Component(:R)
    @test r isa Component{1}
    @test r.type == :R
end

@testset "Inserting Nodes" begin
    bg = RCI()

    c, r, J0 = bg.nodes[[1, 2, 5]]

    bondc0 = getbonds(bg, c, J0)[1]
    bondr0 = getbonds(bg, r, J0)[1]

    tf = Component(:TF, numports = 2)
    insert_node!(bg, bondc0, tf)
    insert_node!(bg, bondr0, EqualFlow())

    @test tf in bg.nodes
    @test nv(bg) == 7
    @test ne(bg) == 6
end

@testset "Merging components" begin
    bg = RCI()

    newC = Component(:C, :newC)
    newR = Component(:R, :newR)
    add_node!(bg, [newC, newR])
    connect!(bg, newC, newR)

    C = getnodes(bg, "C")[1]
    merge_nodes!(bg, C, newC)

    R = getnodes(bg, "R")[1]
    merge_nodes!(bg, R, newR; junction = EqualFlow())

    @test isempty(getnodes(bg, "newC"))
    @test length(getnodes(bg, "1")) == 1
    @test length(getnodes(bg, "0")) == 2
    @test nv(bg) == 7
    @test ne(bg) == 7
end

@testset "Simplifying Junctions" begin
    bg = RCI()
    C, R, I, SS, J0 = bg.nodes

    J0_new_1 = EqualEffort(; name = :new0_1)
    J0_new_2 = EqualEffort(; name = :new0_2)
    insert_node!(bg, (C, J0), J0_new_1)
    insert_node!(bg, (R, J0), J0_new_2)
    connect!(bg, J0_new_1, J0_new_2)

    J1_new_1 = EqualFlow(; name = :new1_1)
    J1_new_2 = EqualFlow(; name = :new1_2)
    add_node!(bg, J1_new_1)
    connect!(bg, J0, J1_new_1)
    insert_node!(bg, (SS, J0), J1_new_2)

    # Removing junction redundancies
    @test length(getnodes(bg, EqualFlow)) == 2
    simplify_junctions!(bg, squash_identical = false)
    @test length(getnodes(bg, EqualFlow)) == 0
    @test nv(bg) == 7
    @test ne(bg) == 7

    # Squashing junction duplicates into a single junction
    simplify_junctions!(bg)
    @test length(getnodes(bg, EqualEffort)) == 1
    @test nv(bg) == 5
    @test ne(bg) == 4
end

@testset "BondGraphNodes" begin
    C = Component(:C)
    bg1 = BondGraph(:BG1)
    bg2 = BondGraph(:BG2)
    bg3 = BondGraph(:BG3)
    main = BondGraph(:Main)

    bgn1 = BondGraphNode(bg1)
    bgn2 = BondGraphNode(bg2)
    bgn3 = BondGraphNode(bg3)

    @test bgn1.bondgraph === bg1
    @test bgn1.type === :BG
    @test bgn1.name === bg1.name
    @test freeports(bgn1) == Bool[]

    add_node!(bg1, C)
    add_node!(bg2, bgn1)
    add_node!(bg3, bgn2)
    add_node!(main, bgn3)

    @test main.BG3.BG2.BG1.C === C

    C2 = Component(:C) # Same name
    add_node!(bg1, C2)
    @test main.BG3.BG2.BG1.C == [C, C2]
end

# @testset "Expose component" begin

# end

# @testset "Nested BondGraphs" begin
#     bg_1 = RCI("RCI_1")
#     bg_2 = RCI("RCI_2")
#     bgn_1 = BondGraphNode(bg_1)
#     bgn_2 = BondGraphNode(bg_2)

#     bg_main = BondGraph("BG_Main")

#     #add_node!(bg_main, [bgn_1, bgn_2])
#     #connect!(bg_main, bgn_1, bgn_2)

# end

max_stars_repo_name: mic-pan/BondGraphs

max_issues_repo_name: mic-pan/BondGraphs

max_forks_repo_name: mic-pan/BondGraphs


text:
using TCX

tcx_sample_file="centry_park_run.tcx"
tcx_centrypark_run_file="centry_park_run.tcx"
tcx_treadmill_run_file="treadmill_run.tcx"
gpx_sample_file="shanghai_marathon_2018.gpx"

include("test_tcx_read_file.jl")
include("test_tcx_read_dir.jl")


text:
module ItemsController

using Genie.Requests
using Genie.Renderer
using Genie.Renderer.Html
using Genie.Renderer.Json
using SearchLight
using Items

export items

# --- Views ---

function items(::Val{:view}, ::Val{:GET})
    html(:items, :item_list; items = all(Item))
end

function items(::Val{:view}, ::Val{:POST})
    d = postpayload()
    @show d
    try
        item = Item(a=d[:a], b=d[:b])
        save(item)
    catch
        return html("Bad Request"; status=400)
    end
    redirect(:get_items)
end

function items(::Val{:view}, ::Val{:GET}, id)
    item = findone(Item; id = id)
    if isnothing(item)
        html("Not Found"; status = 404)
    else
        html(:items, :item; item = item, status = 200)
    end
end

function items(::Val{:view}, ::Val{:POST}, id)
    d = postpayload()
    @show d
    item = findone(Item; id = id)
    if isnothing(item)
        return html("Not Found", status = 404)
    else
        if haskey(d, :delete)
            delete(item)
        else
            try
                item.a = d[:a]
                item.b = d[:b]
                save(item)
            catch
                return html("Bad Request"; status=400)
            end
        end
        redirect(:get_items)
    end
end

# --- API ---

function items(::Val{:api}, ::Val{:GET})
    json(all(Item))
end

function items(::Val{:api}, ::Val{:POST})
    d = jsonpayload()
    @show d
    try
        item = Item(a=d["a"], b=d["b"])
        save(item)
        return json("Created"; status=201)
    catch
        return json("Bad Request"; status=400)
    end
end

function items(::Val{:api}, ::Val{:GET}, id)
    item = findone(Item; id = id)
    if isnothing(item)
        return json("Not Found"; status=404)
    else
        return json(item; status=200)
    end
end

function items(::Val{:api}, ::Val{:PUT}, id)
    d = jsonpayload()
    @show d
    item = findone(Item; id = id)
    if isnothing(item)
        return json("Not Found"; status=404)
    else
        try
            item.a = d["a"]
            item.b = d["b"]
            save(item)
            return json(""; status=200)
        catch
            return json("Bad Request"; status=400)
        end
    end
end

function items(::Val{:api}, ::Val{:DELETE}, id)
    item = findone(Item; id = id)
    if isnothing(item)
        return json("Not Found"; status=404)
    else
        delete(item)
        return json("Deleted"; status=200)
    end
end

end

max_stars_repo_name: jaantollander/genie-webapp-db

max_issues_repo_name: jaantollander/genie-webapp-db

max_forks_repo_name: jaantollander/genie-webapp-db


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

mutable struct IoK8sApiCoreV1EndpointSubset <: SwaggerModel
    addresses::Any # spec type: Union{ Nothing, Vector{IoK8sApiCoreV1EndpointAddress} } # spec name: addresses
    notReadyAddresses::Any # spec type: Union{ Nothing, Vector{IoK8sApiCoreV1EndpointAddress} } # spec name: notReadyAddresses
    ports::Any # spec type: Union{ Nothing, Vector{IoK8sApiCoreV1EndpointPort} } # spec name: ports

    function IoK8sApiCoreV1EndpointSubset(;addresses=nothing, notReadyAddresses=nothing, ports=nothing)
        o = new()
        validate_property(IoK8sApiCoreV1EndpointSubset, Symbol("addresses"), addresses)
        setfield!(o, Symbol("addresses"), addresses)
        validate_property(IoK8sApiCoreV1EndpointSubset, Symbol("notReadyAddresses"), notReadyAddresses)
        setfield!(o, Symbol("notReadyAddresses"), notReadyAddresses)
        validate_property(IoK8sApiCoreV1EndpointSubset, Symbol("ports"), ports)
        setfield!(o, Symbol("ports"), ports)
        o
    end
end # type IoK8sApiCoreV1EndpointSubset

const _property_map_IoK8sApiCoreV1EndpointSubset = Dict{Symbol,Symbol}(Symbol("addresses")=>Symbol("addresses"), Symbol("notReadyAddresses")=>Symbol("notReadyAddresses"), Symbol("ports")=>Symbol("ports"))
const _property_types_IoK8sApiCoreV1EndpointSubset = Dict{Symbol,String}(Symbol("addresses")=>"Vector{IoK8sApiCoreV1EndpointAddress}", Symbol("notReadyAddresses")=>"Vector{IoK8sApiCoreV1EndpointAddress}", Symbol("ports")=>"Vector{IoK8sApiCoreV1EndpointPort}")
Base.propertynames(::Type{ IoK8sApiCoreV1EndpointSubset }) = collect(keys(_property_map_IoK8sApiCoreV1EndpointSubset))
Swagger.property_type(::Type{ IoK8sApiCoreV1EndpointSubset }, name::Symbol) = Union{Nothing,eval(Meta.parse(_property_types_IoK8sApiCoreV1EndpointSubset[name]))}
Swagger.field_name(::Type{ IoK8sApiCoreV1EndpointSubset }, property_name::Symbol) =  _property_map_IoK8sApiCoreV1EndpointSubset[property_name]

function check_required(o::IoK8sApiCoreV1EndpointSubset)
    true
end

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


text:
module Checkboxes

Base.Experimental.@optlevel 1

using Redux
using CImGui

# actions
abstract type AbstractCheckboxAction <: AbstractSyncAction end

"""
    Rename(label, new_label)
Change checkbox's label to `new_label`.
Note that renaming may also change the identifier, please refer to `help?> CImGui.PushID`.
"""
struct Rename <: AbstractCheckboxAction
    label::String
    new_label::String
end

"""
    SetTriggeredTo(label, is_triggered)
Update widget's `is_triggered` state.
"""
struct SetTriggeredTo <: AbstractCheckboxAction
    label::String
    is_triggered::Bool
end

"""
    Toggle(label)
Toggle checkbox.
"""
struct Toggle <: AbstractCheckboxAction
    label::String
end

"""
    Check(label)
Set checkbox to be checked.
"""
struct Check <: AbstractCheckboxAction
    label::String
end

"""
    Uncheck(label)
Set checkbox to be unchecked.
"""
struct Uncheck <: AbstractCheckboxAction
    label::String
end

"""
    SetCheckedTo(label, is_check)
Set checkbox state to `is_check`.
"""
struct SetCheckedTo <: AbstractCheckboxAction
    label::String
    is_check::Bool
end

"""
    AddCheckbox(label::AbstractString, is_check = false)
Add a [`Checkbox`](@ref).
"""
struct AddCheckbox <: AbstractCheckboxAction
    label::String
    is_check::Bool
end
AddCheckbox(label::AbstractString) = AddCheckbox(label, false)

"""
    DeleteCheckbox(label::AbstractString)
Delete the [`Checkbox`](@ref).
"""
struct DeleteCheckbox <: AbstractCheckboxAction
    label::String
end

# state
abstract type AbstractCheckboxState <: AbstractImmutableState end

"""
    Checkboxes.State(label::AbstractString)
    Checkboxes.State(label::AbstractString, is_check)
    Checkboxes.State(label::AbstractString, is_check, is_triggered)
A simple checkbox state which contains a label a.k.a the identifier, a flag value that
marks whether the widget is checked and a flag value `is_triggered` that records the
state of the latest poll events.
"""
struct State <: AbstractCheckboxState
    label::String
    is_check::Bool
    is_triggered::Bool
end
State(label::AbstractString, is_check = false) = State(label, is_check, false)

# reducers
# reducer(state::AbstractState, action::AbstractAction) = state
# reducer(state::Vector{<:AbstractState}, action::AbstractAction) = state
# reducer(state::Dict{String,<:AbstractState}, action::AbstractAction) = state

reducer(s::State, a::Rename) = State(a.new_label, s.is_check, s.is_triggered)
reducer(s::State, a::Toggle) = State(s.label, !s.is_check, s.is_triggered)
reducer(s::State, a::Check) = State(s.label, true, s.is_triggered)
reducer(s::State, a::Uncheck) = State(s.label, false, s.is_triggered)
reducer(s::State, a::SetCheckedTo) = State(s.label, a.is_check, s.is_triggered)
reducer(s::State, a::SetTriggeredTo) = State(s.label, s.is_check, a.is_triggered)

reducer(s::Dict{String,<:AbstractCheckboxState}, a::AbstractCheckboxAction) =
    Dict(k => (get_label(v) == a.label ? reducer(v, a) : v) for (k, v) in s)

reducer(s::Vector{State}, a::AbstractCheckboxAction) = map(s) do s
    get_label(s) == a.label ? reducer(s, a) : s
end

reducer(s::Vector{State}, a::AddCheckbox) =
    State[s..., State(a.label, a.is_check, false)]
reducer(s::Vector{State}, a::DeleteCheckbox) = filter(s -> s.label !== a.label, s)

# helper
"""
    Checkbox(store::AbstractStore, get_state=Redux.get_state) -> Bool
Return `true` when triggered.
`get_state` is a router function that tells how to find the target state from `store`. `chain_action` is for chaining upstream actions.
"""
function Checkbox(store::AbstractStore, get_state=Redux.get_state)
    s = get_state(store)
    check_ref = Ref(s.is_check)
    is_triggered = CImGui.Checkbox(s.label, check_ref)
    dispatch!(store, SetTriggeredTo(s.label, is_triggered))
    is_triggered && dispatch!(store, SetCheckedTo(s.label, check_ref[]))
    return is_triggered
end

get_label(s) = "__REDUX_IMGUI_RESERVED_DUMMY_LABEL"
get_label(s::State) = s.label

is_checked(s::State) = s.is_checked
is_triggered(s::State) = s.is_triggered

end # module


text:
# This z resolution is get from
zs_NCAR_LENS = - [
    0, 1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000,
    11000, 12000, 13000, 14000, 15000, 16000, 17019.681640625,
    18076.12890625, 19182.125, 20349.931640625, 21592.34375, 22923.3125,
    24358.453125, 25915.580078125, 27615.259765625, 29481.470703125,
    31542.373046875, 33831.2265625, 36387.47265625, 39258.046875,
    42498.88671875, 46176.65625, 50370.6875, 55174.91015625, 60699.66796875,
    67072.859375, 74439.8046875, 82960.6953125, 92804.3515625,
    104136.8203125, 117104.015625, 131809.359375, 148290.078125,
    166499.203125, 186301.4375, 207487.390625, 229803.90625, 252990.40625,
    276809.84375, 301067.0625, 325613.84375, 350344.875, 375189.1875,
    400101.15625, 425052.46875, 450026.0625, 475012, 500004.6875,
    525000.9375, 549999.0625
] / 100.0;

#############################################################################
# 2019/05/09
# I want to double the resultion below 2000 meters

# zs_NCAR_LENS[45] = 1863.0
# zs_NCAR_LENS[46] = 2074.9
cut = 45

zs_SSM = []
for i = 1:cut
    append!(zs_SSM, zs_NCAR_LENS[i])
    append!(zs_SSM, (zs_NCAR_LENS[i] + zs_NCAR_LENS[i+1]) / 2.0)
end

for i = cut+1:length(zs_NCAR_LENS)
    append!(zs_SSM, zs_NCAR_LENS[i])
end

#############################################################################

# 2019/05/10
# Make 40 layers (41 points). 5 meters increament in top 100m (20 layers),
# and a linearly increasing increament in the next 20 layers until 1000m.
zs_SSM = zeros(Float64, 41)
zs_SSM[1] = 0.0

dh = zeros(Float64, 40)
dh[ 1:20] .= 5.0
dh[21:40] = 5.0 .+ collect(1:20) * 80.0/21.0
for i = 1:length(dh)
    zs_SSM[i+1] = zs_SSM[i] - dh[i]
end

#############################################################################

# 2019/07/11
# Make 40 layers (41 points). Simply use the same zs as LENS (top 41 pts).

zs_SSM[:] = zs_NCAR_LENS[1:41]

#############################################################################

# 2019/11/03
zs_SSM = collect(Float64, 0:-50:-1000)

#############################################################################

zs_mid_NCAR_LENS = (zs_NCAR_LENS[2:end] + zs_NCAR_LENS[1:end-1] ) / 2.0
zs_mid_SSM = (zs_SSM[2:end] + zs_SSM[1:end-1] ) / 2.0

max_stars_repo_name: meteorologytoday/CESM-diagnostic

max_issues_repo_name: meteorologytoday/CESM-diagnostic

max_forks_repo_name: meteorologytoday/CESM-diagnostic


text:
module EquationOfStateRecipes

using EquationsOfStateOfSolids:
    EquationOfStateOfSolids, EnergyEquation, PressureEquation, BulkModulusEquation
using RecipesBase: @recipe
using Unitful: AbstractQuantity, ustrip, @u_str

import RecipesBase

@recipe function f(
    eos::EquationOfStateOfSolids,
    volumes::AbstractVector = (0.5:0.01:1.1) * eos.param.v0,
)
    if eltype(volumes) <: AbstractQuantity
        u = _yunit(typeof(eos))
        xguide --> "volume (angstrom^3)"
        yguide --> _yprefix(typeof(eos)) * " ($u)"
        x = map(Base.Fix1(ustrip, u"angstrom^3"), volumes)
        y = map(Base.Fix1(ustrip, u) ∘ eos, volumes)
        return x, y
    else
        xguide --> "volume"
        yguide --> _yprefix(typeof(eos))
        return volumes, map(eos, volumes)
    end
end

_yprefix(::Type{<:EnergyEquation}) = "energy"
_yprefix(::Type{<:PressureEquation}) = "pressure"
_yprefix(::Type{<:BulkModulusEquation}) = "bulk modulus"

_yunit(::Type{<:EnergyEquation}) = u"eV"
_yunit(::Type{<:PressureEquation}) = u"GPa"
_yunit(::Type{<:BulkModulusEquation}) = u"GPa"

end


text:
using BenchmarkTools, AMDGPU

DAT = Float64

if DAT==Float64
    DAT_Int = Int64
    sc = 1
elseif DAT==Float32
    DAT_Int = Int32
    sc = 2
end

@inbounds function memcopy_triad!(A, B, C, s)
    ix = (workgroupIdx().x - 1) * workgroupDim().x + workitemIdx().x
    iy = (workgroupIdx().y - 1) * workgroupDim().y + workitemIdx().y
    A[ix,iy] = B[ix,iy] + s*C[ix,iy]
    return
end

@inbounds function memcopy_triad_pow_int!(A, B, C, s, pow_int)
    ix = (workgroupIdx().x - 1) * workgroupDim().x + workitemIdx().x
    iy = (workgroupIdx().y - 1) * workgroupDim().y + workitemIdx().y
    A[ix,iy] = B[ix,iy] + s*C[ix,iy]^pow_int
    return
end

@inbounds function memcopy_triad_pow_float!(A, B, C, s, pow_float)
    ix = (workgroupIdx().x - 1) * workgroupDim().x + workitemIdx().x
    iy = (workgroupIdx().y - 1) * workgroupDim().y + workitemIdx().y
    A[ix,iy] = B[ix,iy] + s*C[ix,iy]^pow_float
    return
end

function diff2D_step!(T2, T, Ci, lam, dt, _dx, _dy)
    ix = (workgroupIdx().x - 1) * workgroupDim().x + workitemIdx().x
    iy = (workgroupIdx().y - 1) * workgroupDim().y + workitemIdx().y
    if (ix>1 && ix<size(T2,1) && iy>1 && iy<size(T2,2))
        @inbounds T2[ix,iy] = T[ix,iy] + dt*(Ci[ix,iy]*(
                              - ((-lam*(T[ix+1,iy] - T[ix,iy])*_dx) - (-lam*(T[ix,iy] - T[ix-1,iy])*_dx))*_dx
                              - ((-lam*(T[ix,iy+1] - T[ix,iy])*_dy) - (-lam*(T[ix,iy] - T[ix,iy-1])*_dy))*_dy ))
    end
    return
end

function run_bench()
    fact      = 24
    nx, ny    = sc*fact*1024, fact*1024
    threads   = (32, 8)
    grid      = (nx, ny)
    A         = ROCArray(zeros(DAT, nx, ny))
    B         = ROCArray( rand(DAT, nx, ny))
    C         = ROCArray( ones(DAT, nx, ny))
    pow_int   = DAT_Int(3)
    pow_float = DAT(3.75)
    s         = rand(DAT)
    lam       = rand(DAT)
    _dx, _dy  = DAT(1.0), DAT(1.0)
    dt        = DAT(1.0/10.0/4.1)
    println("nx, ny, DAT = $(nx), $(ny), $(DAT)")
    # run test 1
    t_it = @belapsed begin wait( @roc groupsize=$threads gridsize=$grid memcopy_triad!($A, $B, $C, $s) ) end
    T_tot = 3*1/1e9*nx*ny*sizeof(DAT)/t_it
    println("T_tot triad2D           = $(round(T_tot,sigdigits=7)) GB/s")
    # run test 2
    t_it = @belapsed begin wait( @roc groupsize=$threads gridsize=$grid memcopy_triad_pow_int!($A, $B, $C, $s, $pow_int) ) end
    T_tot = 3*1/1e9*nx*ny*sizeof(DAT)/t_it
    println("T_tot triad2D pow_int   = $(round(T_tot,sigdigits=7)) GB/s")
    # run test 3
    t_it = @belapsed begin wait( @roc groupsize=$threads gridsize=$grid memcopy_triad_pow_float!($A, $B, $C, $s, $pow_float) ) end
    T_tot = 3*1/1e9*nx*ny*sizeof(DAT)/t_it
    println("T_tot triad2D pow_float = $(round(T_tot,sigdigits=7)) GB/s")
    # run test 4
    t_it = @belapsed begin wait( @roc groupsize=$threads gridsize=$grid diff2D_step!($A, $B, $C, $lam, $dt, $_dx, $_dy) ) end
    T_tot = 3*1/1e9*nx*ny*sizeof(DAT)/t_it
    println("T_tot diffusion 2D      = $(round(T_tot,sigdigits=7)) GB/s")
    return
end

run_bench()

max_stars_repo_name: luraess/JuliaGPUPerf

max_issues_repo_name: luraess/JuliaGPUPerf

max_forks_repo_name: luraess/JuliaGPUPerf


text:
@testset "Eichler" begin
  K, a = Hecke.rationals_as_number_field()
  A = Hecke.AlgQuat(K, K(1), K(1))
  A = Hecke.AlgAss(K, A.mult_table)
  OA = maximal_order(A)
  #b = @inferred Hecke.principal_generator_eichler(2 * OA)
end


text:
"""
    gpuinfo(deviceid::Integer)

Print out detailed information about the GPU with the given `deviceid`.

Heavily inspired by the CUDA sample "deviceQueryDrv.cpp".
"""
function gpuinfo(deviceid::Integer; io::IO=stdout)
    0 <= deviceid <= ngpus() - 1 || throw(ArgumentError("Invalid device id."))
    return gpuinfo(CuDevice(deviceid); io)
end
function gpuinfo(dev::CuDevice=CUDA.device(); io::IO=stdout)
    # query
    mp = nmultiprocessors(dev)
    cores = ncudacores(dev)
    max_clock_rate = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_CLOCK_RATE) ÷ 1000
    mem_clock_rate = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MEMORY_CLOCK_RATE) ÷ 1000
    mem_bus_width = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_GLOBAL_MEMORY_BUS_WIDTH)
    l2cachesize = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_L2_CACHE_SIZE)
    maxTex1D = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE1D_WIDTH)
    maxTex2D_width = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE2D_WIDTH)
    maxTex2D_height = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE2D_HEIGHT)
    maxTex3D_width = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE3D_WIDTH)
    maxTex3D_height = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE3D_HEIGHT)
    maxTex3D_depth = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE3D_DEPTH)
    maxTex1DLayered_width = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE1D_LAYERED_WIDTH
    )
    maxTex1DLayered_layers = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE1D_LAYERED_LAYERS
    )
    maxTex2DLayered_width = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE2D_LAYERED_WIDTH
    )
    maxTex2DLayered_height = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE2D_LAYERED_HEIGHT
    )
    maxTex2DLayered_layers = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAXIMUM_TEXTURE2D_LAYERED_LAYERS
    )
    total_constant_mem = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_TOTAL_CONSTANT_MEMORY)
    shared_mem_per_block = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_SHARED_MEMORY_PER_BLOCK
    )
    regs_per_block = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_REGISTERS_PER_BLOCK)
    warpsize = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_WARP_SIZE)
    max_threads_per_mp = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_THREADS_PER_MULTIPROCESSOR
    )
    max_threads_per_block = CUDA.attribute(
        dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_THREADS_PER_BLOCK
    )
    blockdim_x = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_BLOCK_DIM_X)
    blockdim_y = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_BLOCK_DIM_Y)
    blockdim_z = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_BLOCK_DIM_Z)
    griddim_x = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_GRID_DIM_X)
    griddim_y = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_GRID_DIM_Y)
    griddim_z = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_GRID_DIM_Z)
    texture_align = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_TEXTURE_ALIGNMENT)
    max_mem_pitch = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MAX_PITCH)
    async_engine_count = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_ASYNC_ENGINE_COUNT)
    gpu_overlap = Bool(CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_GPU_OVERLAP))
    kernel_exec_timeout_enabled = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_KERNEL_EXEC_TIMEOUT)
    )
    integrated = Bool(CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_INTEGRATED))
    can_map_host_mem = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_CAN_MAP_HOST_MEMORY)
    )
    concurrent_kernels = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_CONCURRENT_KERNELS)
    )
    surface_alignment = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_SURFACE_ALIGNMENT) > 0
    ecc_enabled = Bool(CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_ECC_ENABLED))
    unified_addressing = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_UNIFIED_ADDRESSING)
    )
    managed_memory = Bool(CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MANAGED_MEMORY))
    compute_preemption = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_COMPUTE_PREEMPTION_SUPPORTED)
    )
    cooperative_launch = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_COOPERATIVE_LAUNCH)
    )
    cooperative_multi_dev_launch = Bool(
        CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_COOPERATIVE_MULTI_DEVICE_LAUNCH)
    )
    pci_domainid = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_PCI_DOMAIN_ID)
    pci_busid = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_PCI_BUS_ID)
    pci_deviceid = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_PCI_DEVICE_ID)
    compute_mode = CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_COMPUTE_MODE)
    comp_modes = [
        "Default (multiple host threads can use ::cudaSetDevice() with device simultaneously)",
        "Exclusive (only one host thread in one process is able to use ::cudaSetDevice() with this device)",
        "Prohibited (no host thread can use ::cudaSetDevice() with this device)",
        "Exclusive Process (many threads in one process is able to use ::cudaSetDevice() with this device)",
        "Unknown",
    ]

    # printing
    println(io, "Device: ", name(dev), " ($dev)")
    println(io, "Total amount of global memory: ", Base.format_bytes(Int(CUDA.totalmem(dev))))
    println(io, "Number of CUDA cores: ", cores)
    println(io, "Number of multiprocessors: ", mp, " ($(cores ÷ mp) CUDA cores each)")
    println(io, "GPU max. clock rate: ", max_clock_rate, " MHz")
    println(io, "Memory clock rate: ", mem_clock_rate, " MHz")
    println(io, "Memory bus width: ", mem_bus_width, "-bit")
    println(io, "L2 cache size: ", Base.format_bytes(l2cachesize))
    println(io, "Max. texture dimension sizes (1D): $maxTex1D")
    println(io, "Max. texture dimension sizes (2D): $maxTex2D_width, $maxTex2D_height")
    println(io,
        "Max. texture dimension sizes (3D): $maxTex3D_width, $maxTex3D_height, $maxTex3D_depth",
    )
    println(io,
        "Max. layered 1D texture size: $(maxTex1DLayered_width) ($(maxTex1DLayered_layers) layers)",
    )
    println(io,
        "Max. layered 2D texture size: $(maxTex2DLayered_width), $(maxTex2DLayered_height) ($(maxTex2DLayered_layers) layers)",
    )
    println(io, "Total amount of constant memory: ", Base.format_bytes(total_constant_mem))
    println(io,
        "Total amount of shared memory per block: ", Base.format_bytes(shared_mem_per_block)
    )
    println(io, "Total number of registers available per block: ", regs_per_block)
    println(io, "Warp size: ", warpsize)
    println(io, "Max. number of threads per multiprocessor: ", max_threads_per_mp)
    println(io, "Max. number of threads per block: ", max_threads_per_block)
    println(io,
        "Max. dimension size of a thread block (x,y,z): $(blockdim_x), $(blockdim_y), $(blockdim_z)",
    )
    println(io,
        "Max. dimension size of a grid size (x,y,z): $(griddim_x), $(griddim_y), $(griddim_z)",
    )
    println(io, "Texture alignment: ", Base.format_bytes(texture_align))
    println(io, "Maximum memory pitch: ", Base.format_bytes(max_mem_pitch))
    println(io,
        "Concurrent copy and kernel execution: ",
        gpu_overlap ? "Yes" : "No",
        " with $(async_engine_count) copy engine(s)",
    )
    println(io, "Run time limit on kernels: ", kernel_exec_timeout_enabled ? "Yes" : "No")
    println(io, "Integrated GPU sharing host memory: ", integrated ? "Yes" : "No")
    println(io, "Support host page-locked memory mapping: ", can_map_host_mem ? "Yes" : "No")
    println(io, "Concurrent kernel execution: ", concurrent_kernels ? "Yes" : "No")
    println(io, "Alignment requirement for surfaces: ", surface_alignment ? "Yes" : "No")
    println(io, "Device has ECC support: ", ecc_enabled ? "Yes" : "No")
    println(io, "Device supports Unified Addressing (UVA): ", unified_addressing ? "Yes" : "No")
    println(io, "Device supports managed memory: ", managed_memory ? "Yes" : "No")
    println(io, "Device supports compute preemption: ", compute_preemption ? "Yes" : "No")
    println(io, "Supports cooperative kernel launch: ", cooperative_launch ? "Yes" : "No")
    println(io,
        "Supports multi-device co-op kernel launch: ",
        cooperative_multi_dev_launch ? "Yes" : "No",
    )
    println(io,
        "Device PCI domain ID / bus ID / device ID: $(pci_domainid) / $(pci_busid) / $(pci_deviceid)",
    )
    println(io, "Compute mode: ", comp_modes[compute_mode + 1])

    return nothing
end

"""
Query peer-to-peer (i.e. inter-GPU) access support.
"""
function gpuinfo_p2p_access(; io::IO=stdout)
    # check p2p access
    ndevs = ngpus()
    if ndevs <= 1
        println(io,"Only a single GPU available.")
    else
        mat_p2p_access_supported = Matrix{Bool}(undef, ndevs, ndevs)
        mat_p2p_can_access = Matrix{Bool}(undef, ndevs, ndevs)
        mat_p2p_atomic_supported = Matrix{Bool}(undef, ndevs, ndevs)
        for i in 1:ndevs
            dev_i = CuDevice(i - 1)
            for j in 1:ndevs
                dev_j = CuDevice(j - 1)
                if i != j
                    p2p_access_supported = Bool(
                        CUDA.p2p_attribute(
                            dev_i, dev_j, CUDA.CU_DEVICE_P2P_ATTRIBUTE_ACCESS_SUPPORTED
                        ),
                    )
                    p2p_can_access = Bool(CUDA.can_access_peer(dev_i, dev_j))
                    p2p_atomic_supported = Bool(
                        CUDA.p2p_attribute(
                            dev_i,
                            dev_j,
                            CUDA.CU_DEVICE_P2P_ATTRIBUTE_NATIVE_ATOMIC_SUPPORTED,
                        ),
                    )
                    mat_p2p_atomic_supported[i, j] = p2p_atomic_supported
                    mat_p2p_access_supported[i, j] = p2p_access_supported
                    mat_p2p_can_access[i, j] = p2p_can_access
                    # p2p_performance_rank = CUDA.p2p_attribute(dev_i, dev_j, CUDA.CU_DEVICE_P2P_ATTRIBUTE_PERFORMANCE_RANK)
                else
                    mat_p2p_atomic_supported[i, i] = false
                    mat_p2p_access_supported[i, i] = false
                    mat_p2p_can_access[i, j] = false
                end
            end
        end

        printstyled(io,"P2P Access Supported:\n"; bold=true)
        show(io, "text/plain",mat_p2p_access_supported)
        println(io)
        println(io)
        if mat_p2p_access_supported != mat_p2p_can_access
            printstyled(io,"P2P Can Access:\n"; bold=true)
            show(io, "text/plain",mat_p2p_can_access)
            println(io)
            println(io)
        end
        printstyled(io,"P2P Atomic Supported:\n"; bold=true)
        show(io, "text/plain",mat_p2p_atomic_supported)
        println(io)
        println(io)
    end
    return nothing
end

ngpus() = length(CUDA.devices())
function nmultiprocessors(dev::CuDevice=CUDA.device())
    return CUDA.attribute(dev, CUDA.CU_DEVICE_ATTRIBUTE_MULTIPROCESSOR_COUNT)
end
ncudacores(deviceid::Integer) = ncudacores(CuDevice(deviceid))
function ncudacores(device::CuDevice=CUDA.device())
    capver = CUDA.capability(device)
    return ncudacores(capver.major, capver.minor, nmultiprocessors(device))
end
function ncudacores(major, minor, mp)
    # based on https://stackoverflow.com/questions/32530604/how-can-i-get-number-of-cores-in-cuda-device
    # helper_cuda_drvapi provides something like https://github.com/LinkedInAttic/datacl/blob/master/approxalgos/GPU_Work_Final2/bussAnal/filter/lib/helper_cuda_drvapi.h#L82 but is header only
    cores = 0
    if major == 2 # Fermi
        if minor == 1
            cores = mp * 48
        else
            cores = mp * 32
        end
    elseif major == 3 # Kepler
        cores = mp * 192
    elseif major == 5 # Maxwell
        cores = mp * 128
    elseif major == 6 # Pascal
        if (minor == 1) || (minor == 2)
            cores = mp * 128
        elseif minor == 0
            cores = mp * 64
        else
            error("Unknown device type")
        end
    elseif major == 7 # Volta and Turing
        if (minor == 0) || (minor == 5)
            cores = mp * 64
        else
            error("Unknown device type")
        end
    elseif major == 8 # Ampere
        if minor == 0
            cores = mp * 64
        elseif minor == 6
            cores = mp * 128
        else
            error("Unknown device type")
        end
    else
        error("Unknown device type")
    end
    return cores
end

function ntensorcores(device::CuDevice=CUDA.device())
    capver = CUDA.capability(device)
    return ntensorcores(capver.major, capver.minor, nmultiprocessors(device))
end
function ntensorcores(major, minor, mp)
    # based on https://en.wikipedia.org/wiki/CUDA
    if major == 8
        return 4 * mp
    elseif major == 7
        return 8 * mp
    else
        return 0
    end
end

"""
List the available GPUs.
"""
function gpus(; io::IO=stdout)
    # Based on https://github.com/JuliaGPU/CUDA.jl/blob/ca77d1828f3bc0df34501de848c7a13f1df0b1fe/src/utilities.jl#L69
    devs = devices()
    if isempty(devs)
        println(io,"No CUDA-capable devices.")
    elseif length(devs) == 1
        println(io,"1 device:")
    else
        println(io,length(devs), " devices:")
    end
    for (i, dev) in enumerate(devs)
        if has_nvml()
            mig = uuid(dev) != parent_uuid(dev)
            nvml_gpu = NVML.Device(parent_uuid(dev))
            nvml_dev = NVML.Device(uuid(dev); mig)

            str = NVML.name(nvml_dev)
            cap = NVML.compute_capability(nvml_gpu)
            mem = NVML.memory_info(nvml_dev)
        else
            str = name(dev)
            cap = capability(dev)
            mem = device!(dev) do
                # this requires a device context, so we prefer NVML
                (free=available_memory(), total=total_memory())
            end
        end
        println(io,
            "  $(i-1): $str (sm_$(cap.major)$(cap.minor), $(Base.format_bytes(mem.free)) / $(Base.format_bytes(mem.total)) available)",
        )
    end
end


text:
using DashBootstrapComponents, DashHtmlComponents

progress = html_div([
    dbc_progress(value = 80, id = "animated-progress", animated = false, striped = true),
    dbc_button(
        "Toggle animation",
        id = "animation-toggle",
        className = "mt-3",
        n_clicks = 0,
    ),
]);

callback!(
    app,
    Output("animated-progress", "animated"),
    Input("animation-toggle", "n_clicks"),
    State("animated-progress", "animated"),
) do n, animated
    return n > 0 ? animated == 0 : animated
end;

max_stars_repo_name: glsdown/dash-bootstrap-components

max_issues_repo_name: glsdown/dash-bootstrap-components

max_forks_repo_name: glsdown/dash-bootstrap-components


text:
using LazySets, SparseArrays

import LazySets: genmat, dim, linear_map, minkowski_sum, cartesian_product, quadratic_map

"""
    SparsePolynomialZonotope

    Type that represents a sparse polynomial zonotope.

Given a constant offset ``c ∈ \\mathbb{R}^n``, a generator matrix ``G ∈ \\mathbb{R}^{n \times h}``
and an exponent matrix ``\\mathbb{N}^{p×h}_{≥0}``, a sparse polynomial zonotope
``\\mathcal{PZ} ⊂ \\mathbb{R}^n``

### Fields

- `c` -- constant offset vector
- `G` -- generators matrix
- `E` -- exponents matrix

### Notes

Sparse polynomial zonotopes were introduced in N. Kochdumper and M. Althoff.
*Sparse Polynomial Zonotopes: A Novel Set Representation for Reachability Analysis*.
Transactions on Automatic Control, 2021.
"""
struct SparsePolynomialZonotope{N, VN<:AbstractVector{N}, MNG<:AbstractMatrix{N}, MNE<:AbstractMatrix{N}} <: LazySet{N}
    c::VN
    G::MNG
    E::MNE

    function SparsePolynomialZonotope(c::VN, G::MNG, E::MNE) where {N, VN<:AbstractVector{N}, MNG<:AbstractMatrix{N}, MNE<:AbstractMatrix{N}}
        @assert length(c) == size(G, 1) "the length of the offset vector, $(length(c)), is not compatible with a generators matrix of size $(size(G))"
        @assert size(G, 2) == size(E, 2) "a generators matrix of size $(size(G)) is not compatible with a exponent matrix of size $(size(E))"

        return new{N, VN, MNG, MNE}(c, G, E)
    end
end

offset(PZ::SparsePolynomialZonotope) = PZ.c
genmat(PZ::SparsePolynomialZonotope) = PZ.G

dim(PZ::SparsePolynomialZonotope) = size(PZ.G, 1)

function linear_map(M::AbstractMatrix, SP::SparsePolynomialZonotope)
    return SparsePolynomialZonotope(M * SP.c, M * SP.G, E)
end

function minkowski_sum(SP1::SparsePolynomialZonotope, SP2::SparsePolynomialZonotope)
    c = SP1.c + SP2.c
    G = hcat(SP1.G, SP2.G)
    E = _blockdiag(SP1.E, SP2.E)
    return SparsePolynomialZonotope(c, G, E)
end

_blockdiag(A::SparseMatrixCSC, B::SparseMatrixCSC) = blockdiag(A, B)
_blockdiag(A::AbstractMatrix, B::AbstractMatrix) = [A zeros(size(A, 1), size(B, 2)); zeros(size(B, 1), size(A, 2)) B]

function cartesian_product(SP1::SparsePolynomialZonotope, SP2::SparsePolynomialZonotope)
    c = vcat(SP1.c, SP2.c)
    G = _blockdiag(SP1.G, SP2.G)
    E = _blockdiag(SP1.E, SP2.E)
    return SparsePolynomialZonotope(c, G, E)
end

# TODO
function quadratic_map(Q::Vector{<:AbstractMatrix}, PZ::PolynomialZonotope)
    c̄ = [dot(c, Qi, c) for Qi in Q]
    #Ḡ = vcat([dot(c, (Qi+Qi')*G) for Qi in Q])
    # ...
end


text:
# Implementation of an elementary cellular automata
# according to https://mathworld.wolfram.com/ElementaryCellularAutomaton.html
# uses random initialization or single cell initialization and uses wolfram codes to specify the rule
# synchronous update of the 1D lattice

using Agents, Random
using CairoMakie
using InteractiveDynamics
using CSV

"""
The automaton living in a 1D space
"""
mutable struct Cell <: AbstractAgent
    id::Int
    pos::Dims{1}
    status::Int # either 0 or 1, where 1 is 'alive'
end

"""
Returns an array with the rule for the next status configurations, 1-indexed.

Thus, the first position corresponds to the rule for the cell neighborhood status (0,0,0)
and the last to (1,1,1)
"""
function rule_from_code(wolfram_code)
    return digits(wolfram_code, base=2, pad=8)
end

"""
Takes the status of a neighborhood and returns the corresponding
index in the model rule (0,0,0) -> 1 (0,1,0) -> 3
"""
function configuration_index(cell_statuses)
    # takes the tuples and forms a string (0,1,1) -> "011"
    binary_code = string(cell_statuses...)
    # turns the string binary code into the integer
    index = parse(Int, binary_code, base=2) + 1
    return index
end

"""
Given a cell checks its neighbors and
decides the next status of the cell based on the model rule
"""
function next_status(cell, model)
    neighbors = collect(nearby_agents(cell, model))
    cell_statuses = (c.status for c in [neighbors[1], cell, neighbors[2]])
    index = configuration_index(cell_statuses)
    return model.rule[index]
end

"""
Initializes the ABM

# Arguments
- `n_cells::Int`: total of cell automaton
- `wolfram_code::Int`: wolfram code for the rule
- `random_start::Bool`: random start or single cell at center
- `seed::Int`: random seed
"""
function build_model(; n_cells = 101, wolfram_code=30, random_start=true, seed = 30)
    space = GridSpace((n_cells,); metric=:chebyshev)

    properties = Dict(:rule => rule_from_code(wolfram_code),)
    model = ABM(
        Cell,
        space; properties,
        rng = MersenneTwister(seed))

	if random_start
		for x in 1:n_cells
			cell = Cell(nextid(model), (x,), rand([0,1]))
			add_agent_pos!(cell, model)
		end
	else
		for x in 1:n_cells
			initial_state = x == div(n_cells, 2) + 1 ? 1 : 0
			cell = Cell(nextid(model), (x,), initial_state)
        	add_agent_pos!(cell, model)
		end
	end

    return model
end

"""
Dummy update of the cells
"""
function cell_step!(_, _)

end

"""
Performs a synchronous update of all cells
"""
function ca_step!(model)
    new_statuses = fill(0, nagents(model))
    for agent in allagents(model)
        new_statuses[agent.id] = next_status(agent, model)
    end

    for id in allids(model)
        model[id].status = new_statuses[id]
    end
end

# Initialize model
model = build_model(n_cells=101, wolfram_code=30, random_start=false)

# Runs the model and collects data
data, _ = run!(model, cell_step!, ca_step!, 100; adata=[:status]);

# The data contains the step, id/position and status of the cell (1/0)
data

CSV.write("data/ca_data.csv", data);

# Lets plot the time evolution in the y axis
function plot_automata(data)
    fig, ax = heatmap(data.id, data.step, data.status, colormap=:Blues_3)
    ax.yreversed = true
    return fig
end

plot_automata(data)


text:
using StatsDiscretizations
using Test
using LinearAlgebra
using Infinity
using Infinity.Utils
using SplitApplyCombine

@testset "Distributions Intervals" begin
    include("test_distributions_intervals.jl")
end

@testset "Discretizers" begin
    include("test_discretizers.jl")
end

@testset "Discretized Functions" begin
    include("test_discretized_functions.jl")
end


text:
# helper functions for vector inputs

@vectorize_1arg Date isweekend

function isbday(hc::HolidayCalendar, dt::Vector{Date})
    result = Array(Bool, length(dt))

    for i in eachindex(dt)
        result[i] = isbday(hc, dt[i])
    end

    return result
end

function isbday(hc::Vector{HolidayCalendar}, dt::Vector{Date})
    const l_hc = length(hc)
    const l_dt = length(dt)

    if l_hc != l_dt
        error("Input vectors must have the same size. $(l_hc) != $(l_dt)")
    end

    result = Array(Bool, length(dt))

    for i in 1:l_hc
        result[i] = isbday(hc[i], dt[i])
    end

    return result
end

isbday(calendar, dt::Vector{Date}) = isbday(convert(HolidayCalendar, calendar), dt)
isbday{A<:AbstractArray}(calendars::A, dt::Vector{Date}) = isbday(convert(Vector{HolidayCalendar}, calendars), dt)

function tobday(hc::HolidayCalendar, dt::Vector{Date}; forward::Bool = true)
    result = Array(Date, length(dt))

    for i in eachindex(dt)
        result[i] = tobday(hc, dt[i]; forward=forward)
    end

    return result
end

function tobday(hc::Vector{HolidayCalendar}, dt::Vector{Date}; forward::Bool = true)
    const l_hc = length(hc)
    const l_dt = length(dt)

    if l_hc != l_dt
        error("Input vectors must have the same size. $(l_hc) != $(l_dt)")
    end

    result = Array(Date, l_hc)

    for i in 1:l_hc
        result[i] = tobday(hc[i], dt[i]; forward=forward)
    end

    return result
end

tobday(calendar, dt::Vector{Date}; forward::Bool = true) = tobday(convert(HolidayCalendar, calendar), dt; forward=forward)
tobday{A<:AbstractArray}(calendars::A, dt::Vector{Date}; forward::Bool = true) = tobday(convert(Vector{HolidayCalendar}, calendars), dt; forward=forward)

function bdays(hc::HolidayCalendar, dt0::Vector{Date}, dt1::Vector{Date})
    const l0 = length(dt0)
    const l1 = length(dt1)

    if l0 != l1
        error("Input vectors must have the same size. $(l0) != $(l1)")
    end

    result = Array(Day, l0)

    for i in 1:l0
        result[i] = bdays(hc, dt0[i], dt1[i])
    end

    return result
end

function bdays(hc::Vector{HolidayCalendar}, dt0::Vector{Date}, dt1::Vector{Date})
    const l_hc = length(hc)
    const l0 = length(dt0)
    const l1 = length(dt1)

    if l_hc != l0 || l0 != l1
        error("Input vectors must have the same size. $(l_hc), $(l0), $(l1)")
    end

    result = Array(Day, l0)

    for i in 1:l0
        result[i] = bdays(hc[i], dt0[i], dt1[i])
    end

    return result
end

bdays(calendar, dt0::Vector{Date}, dt1::Vector{Date}) = bdays(convert(HolidayCalendar, calendar), dt0, dt1)
bdays{A<:AbstractArray}(calendars::A, dt0::Vector{Date}, dt1::Vector{Date}) = bdays(convert(Vector{HolidayCalendar}, calendars), dt0, dt1)

function advancebdays(hc::HolidayCalendar, dt::Date, bdays_count_vec::Vector{Int})
    l = length(bdays_count_vec)
    result = Array(Date, l)
    for i in 1:l
        result[i] = advancebdays(hc, dt, bdays_count_vec[i])
    end
    return result
end

advancebdays(calendar, dt::Date, bdays_count_vec::Vector{Int}) = advancebdays(convert(HolidayCalendar, calendar), dt, bdays_count_vec)
advancebdays(hc, dt::Date, bdays_range::UnitRange{Int}) = advancebdays(hc, dt, collect(bdays_range))


text:
using LinearOptimizationSolver
using Test
using Random
using LinearAlgebra

include("test_simplex.jl")
include("test_canonical.jl")
include("test_interior_points.jl")
include("test_crossover.jl")


text:
@testset "Test Function Spaces" begin

using jFEMTools
import Tensors: Vec, gradient, dot, norm, Tensor
using LinearAlgebra
jF = jFEMTools

# Load mesh
mesh = unitSquareMesh2(TriangleCell, (2,2));

# ### Define Function Spaces
dim = 2
finiteElement = ContinuousLagrange(:Triangle,1)
Wh = FEMFunctionSpace(mesh, finiteElement, 1)    #Scalar Space
Vh = FEMFunctionSpace(mesh, finiteElement, 2)    #Vector Space

# Basic Test
@test jF.getnlocaldofs(Vh) == 6
@test jF.getnlocaldofs(Wh) == 3
@inbounds for (cellcount, cell) in enumerate(jF.CellIterator(mesh))
    jF.reinit!(Wh, cell)
    jF.reinit!(Vh, cell)
end

# Wh Scalar Space testset
@test jF.getnquadpoints(Wh) == 3
@test jF.getnbasefunctions(Wh) == 3
@test jF.getdetJdV(Wh, 1) ≈ 1/24
@test jF.geometric_value(Wh, 1, 1) ≈ 2/3
@test jF.getdim(Wh) == 2
@test jF.shape_value(Wh,1,1) ≈ 2/3
@test jF.shape_gradient(Wh,1,1) ≈ Vec{2}((0.0,-2.0))
@test jF.shape_divergence(Wh,1,1) ≈ -2.0

# Vh Vector Space testset
@test jF.getnquadpoints(Vh) == 3
@test jF.getnbasefunctions(Vh) == 6
@test jF.getdetJdV(Vh, 1) ≈ 1/24
@test jF.geometric_value(Vh, 1, 1) ≈ 2/3
@test jF.getdim(Vh) == 2
@test jF.shape_value(Vh,1,1) ≈ [2/3, 0.0]
@test jF.shape_gradient(Vh,1,1) ≈ Tensor{2,2}([0.0 -2.0;0.0 0.0])
@test jF.shape_divergence(Vh,1,2) ≈ 2.0

end

max_stars_repo_name: Paulms/jFEMToools

max_issues_repo_name: Paulms/jFEMToools

max_forks_repo_name: Paulms/jFEMToools