text:
function _print_result_value(f, si::SamplesInfo; indent=0, kwargs...)
    print(f, "$(repeat(' ', indent))SamplesInfo(;")
    print(f, join(["$k=$(repr(v))" for (k, v) in pairs(si)], ", "))
    print(f, ")")
end

function _print_result_value(f, p::Pair{A,Vector{B}}; indent=0, tab=4) where {A,B}
    _print_result_value(f, first(p); indent=indent, tab=tab)
    println(f, " => [")
    for x in last(p)
        _print_result_value(f, x; indent=(indent + tab), tab=tab)
        println(f, ",")
    end
    print(f, "$(repeat(' ', indent))]")
end

function _print_result_value(f, v::Vector; indent=0, tab=4)
    println(f, "$(repeat(' ', indent))[")
    for x in v
        _print_result_value(f, x; indent=indent + tab, tab=tab)
        println(f, ",")
    end
    print(f, "$(repeat(' ', indent))]")
end

_print_result_value(f, e; indent=0, kwargs...) = print(f, "$(repeat(' ', indent))$(repr(e))")

function print_result(f, r, label; indent=0, tab=4)
    println(f, "$(repeat(' ', indent))($(repeat(' ', indent-1))# $label")
    for (field, value) in pairs(r)
        println(f, "$(repeat(' ', indent + tab))$field=")
        _print_result_value(f, value; indent=2 * tab + indent)
        println(f, ",")
    end
    print(f, "$(repeat(' ', indent)))")
end

function prettyprint_diagnostic_info(f::IO, results)
    println(f, "[")
    for (i, result) in enumerate(results)
        print_result(f, result, i; indent=4)
        println(f, ",")
    end
    print(f, "]")
end

_edf_headers(r) = vcat(map(last, r.header_map)..., r.unextracted_edf_headers)

_groupby(r) = Set((h.label, h.transducer_type, h.physical_dimension)
                  for h in _edf_headers(r))

"""
    prettyprint_diagnostic_info(filename_base::String, diagnostics; dedup=true)

Write `\$(filename_base).out` with a julia-readable pretty-printing
of an iterable of `edf_to_onda_samples` or `edf_header_to_onda_samples_info`
diagnostic `NamedTuples`.

If `dedup == true` (default), the diagnostics will be de-duplicated
based on the set of original EDF headers they each contain. Each
input result that does not already have a `:nrecordings` field
will get one equal to the number of "duplicate" recordings.

This is used by `test/import.jl`, and an example output can be found at
`test/test_edf_to_samples_info.in`.
"""
function prettyprint_diagnostic_info(filename_base::String, results; dedup=true)
    if dedup
        counts = countmap(map(_groupby, filter(r -> !hasproperty(r, :nrecordings), results)))
        foreach(r -> counts[_groupby(r)] = r.nrecordings, filter(r -> hasproperty(r, :nrecordings), results))
        grouped = Dict(_groupby(r) => r for r in results)
        results = collect((;r..., nrecordings=counts[group]) for (group, r) in pairs(grouped))
    end
    filename = "$(filename_base).out"
    open(filename, "w") do f
        print(f, "$(filename_base) = ")
        prettyprint_diagnostic_info(f, results)
        println(f)
    end
end

mock_edf_signal(header, n_records) = EDF.Signal(header, rand(Int16, n_records * header.samples_per_record))

"""
    mock_edf(diagnostic_info)

Turn a diagnostic info `NamedTuple` as returned by `edf_to_onda_samples` or
`edf_header_to_onda_samples_info` and contained in `test_edf_to_samples_info.in`
into a mock edf, for testing purposes.
"""
function mock_edf(result)
    original_edf_headers = OndaEDF._edf_headers(result)
    n_records = 100
    sample_rate = 256
    samples_per_record = first(original_edf_headers).samples_per_record
    edf_signals = Union{EDF.AnnotationsSignal,EDF.Signal}[mock_edf_signal(t, n_records) for t in original_edf_headers]
    edf_header = EDF.FileHeader("0", "", "", DateTime("2014-10-27T22:24:28"), true, n_records, samples_per_record / sample_rate)
    edf = EDF.File((io = IOBuffer(); close(io); io), edf_header, edf_signals)
    return edf
end


text:
# Test of the transmission and distribution combined model

## Import packages and choose a solver

using Printf
import PowerModels; const _PM = PowerModels
import PowerModelsACDC; const _PMACDC = PowerModelsACDC
import FlexPlan; const _FP = FlexPlan
import Cbc
optimizer = _FP.optimizer_with_attributes(Cbc.Optimizer, "logLevel"=>0)

# Include profile data from FlexPlan repository; you can of course also use your own data
include("../test/io/create_profile.jl")

## Input parameters

# Number of hourly optimization periods
number_of_hours = 12

# PowerModels, PowerModelsACDC and FlexPlan settings
s = Dict("output" => Dict("branch_flows" => true), "conv_losses_mp" => false)

## Transmission network instance (all data preparation except for `make_multinetwork` call)

t_file = "./test/data/case6/case6_2030.m" # Input case for transmission network

t_data = _FP.parse_file(t_file)
_FP.add_dimension!(t_data, :hour, number_of_hours)
_FP.add_dimension!(t_data, :scenario, Dict(1 => Dict{String,Any}("probability"=>1)), metadata = Dict{String,Any}("mc"=>true))
_FP.add_dimension!(t_data, :year, 1; metadata = Dict{String,Any}("scale_factor"=>1))
_FP.add_dimension!(t_data, :sub_nw, 1)
_FP.scale_data!(t_data)
t_data, t_loadprofile, t_genprofile = create_profile_data_italy!(t_data)
t_time_series = create_profile_data(number_of_hours, t_data, t_loadprofile, t_genprofile)

## Distribution network instance 1 (all data preparation except for `make_multinetwork` call)

d_file     = "test/data/cigre_mv_eu/cigre_mv_eu.m" # Input case for distribution networks
scale_load = 1.0 # Scaling factor of loads
scale_gen  = 1.0 # Scaling factor of generators

d_data_1 = _FP.parse_file(d_file)
_FP.add_dimension!(d_data_1, :hour, number_of_hours)
_FP.add_dimension!(d_data_1, :scenario, Dict(1 => Dict{String,Any}("probability"=>1)))
_FP.add_dimension!(d_data_1, :year, 1; metadata = Dict{String,Any}("scale_factor"=>1))
_FP.add_dimension!(d_data_1, :sub_nw, 1)
_FP.shift_ids!(d_data_1, _FP.dim_length(t_data))
_FP.scale_data!(d_data_1)
d_time_series = create_profile_data_cigre(d_data_1, number_of_hours; scale_load, scale_gen)
_FP.add_td_coupling_data!(t_data, d_data_1; t_bus = 1, sub_nw = 1) # The first distribution network is connected to bus 1 of transmission network.

## Distribution network instance 2 (all data preparation except for `make_multinetwork` call)

d_data_2 = deepcopy(d_data_1) # For simplicity, here a second distributon network is generated by duplicating the existing one
_FP.shift_ids!(d_data_2, _FP.dim_length(t_data))
_FP.add_td_coupling_data!(t_data, d_data_2; t_bus = 2, sub_nw = 1) # The second distribution network is connected to bus 2 of transmission network.

## Multinetwork data preparation

t_mn_data = _FP.make_multinetwork(t_data, t_time_series)
d_data_1 = _FP.make_multinetwork(d_data_1, d_time_series)
d_data_2 = _FP.make_multinetwork(d_data_2, d_time_series)
d_mn_data = _FP.merge_multinetworks!(d_data_1, d_data_2, :sub_nw)

## Solve problem

result = _FP.flex_tnep(t_mn_data, d_mn_data, _PM.DCPPowerModel, _FP.BFARadPowerModel, optimizer; setting=s)
@assert result["termination_status"] ∈ (_PM.OPTIMAL, _PM.LOCALLY_SOLVED) "$(result["optimizer"]) termination status: $(result["termination_status"])"

## Analyze result

t_nws = _FP.dim_length(t_mn_data)
d_nws = _FP.dim_length(d_mn_data)
sub_nws = _FP.dim_length(d_mn_data, :sub_nw)

printstyled("\n\n====================   Transmission network, first period   ====================\n\n", bold=true, color=:cyan)
_PM.print_summary(result["solution"]["nw"]["1"])

for sub_nw in 1:sub_nws
    printstyled("\n\n===================   Distribution network $sub_nw, first period   ===================\n\n", bold=true, color=:cyan)
    _PM.print_summary(result["solution"]["nw"]["$(sub_nw*t_nws+1)"])
end

printstyled("\n\n==========================   Power exchange at PCCs   ==========================\n\n", bold=true, color=:cyan)
println("Power in MW and MVar, positive if from transmission to distribution\n")
print("period ")
for sub_nw in 1:sub_nws
    @printf("%13s%10s", "p_dist$sub_nw", "q_dist$sub_nw")
end
println()
for t_nw in 1:t_nws
    @printf("%6i:", t_nw)
    for sub_nw in 1:sub_nws
        d_nw = sub_nw * t_nws + t_nw
        t_gen = _FP.dim_prop(d_mn_data, :sub_nw, sub_nw, "t_gen")
        t_res = result["solution"]["nw"]["$t_nw"]
        d_res = result["solution"]["nw"]["$d_nw"]
        t_p_in = t_res["gen"]["$t_gen"]["pg"] * t_res["baseMVA"]
        d_p_in = d_res["td_coupling"]["p"] * d_res["baseMVA"]
        @assert d_p_in ≈ -t_p_in
        d_q_in = d_res["td_coupling"]["q"] * d_res["baseMVA"]
        @printf("%13.3f%10.3f", d_p_in, d_q_in)
    end
    println()
end

ext: jl


text:
# Standard LP interface
require(joinpath(Pkg.dir("MathProgBase"),"src","MathProgSolverInterface.jl"))
importall MathProgSolverInterface

###############################################################################
# Solver objects
export IpoptSolver
immutable IpoptSolver <: AbstractMathProgSolver
  options
end
IpoptSolver(;kwargs...) = IpoptSolver(kwargs)

type IpoptMathProgModel <: AbstractMathProgModel
  inner::Any
  options
end
function IpoptMathProgModel(;options...)
  return IpoptMathProgModel(nothing,options)
end
model(s::IpoptSolver) = IpoptMathProgModel(;s.options...)
export model

###############################################################################
# Begin interface implementation
function loadproblem!(model::IpoptMathProgModel, A, l, u, c, lb, ub, sense)
  Asparse = convert(SparseMatrixCSC{Float64,Int32}, A)
  n = int(Asparse.n)
  m = int(Asparse.m)
  nnz = int(length(Asparse.rowval))
  c_correct = float(c)
  if sense == :Max
    c_correct .*= -1.0
  end

  # Objective callback
  function eval_f(x)
    return dot(x,c_correct)
  end

  # Objective gradient callback
  function eval_grad_f(x, grad_f)
    for j = 1:n
      grad_f[j] = c_correct[j]
    end
  end

  # Constraint value callback
  function eval_g(x, g)
    g_val = A*x
    for i = 1:m
      g[i] = g_val[i]
    end
  end

  # Jacobian callback
  function eval_jac_g(x, mode, rows, cols, values)
    if mode == :Structure
      # Convert column wise sparse to triple format
      idx = 1
      for col = 1:n
        for pos = Asparse.colptr[col]:(Asparse.colptr[col+1]-1)
          rows[idx] = Asparse.rowval[pos]
          cols[idx] = col
          idx += 1
        end
      end
    else
      # Values
      idx = 1
      for col = 1:n
        for pos = Asparse.colptr[col]:(Asparse.colptr[col+1]-1)
          values[idx] = Asparse.nzval[pos]
          idx += 1
        end
      end
    end
  end

  x_L = float(l)
  x_U = float(u)
  g_L = float(lb)
  g_U = float(ub)
  model.inner = createProblem(n, x_L, x_U, m, g_L, g_U, nnz, 0,
                              eval_f, eval_g, eval_grad_f, eval_jac_g, nothing)
  model.inner.sense = sense
  addOption(model.inner, "jac_c_constant", "yes")
  addOption(model.inner, "jac_d_constant", "yes")
  addOption(model.inner, "hessian_constant", "yes")
  addOption(model.inner, "hessian_approximation", "limited-memory")
  addOption(model.inner, "mehrotra_algorithm", "yes")
  for (name,value) in model.options
    addOption(model.inner, string(name), value)
  end
end

getsense(m::IpoptMathProgModel) = m.inner.sense
numvar(m::IpoptMathProgModel) = m.inner.n
numconstr(m::IpoptMathProgModel) = m.inner.m
optimize!(m::IpoptMathProgModel) = solveProblem(m.inner)
function status(m::IpoptMathProgModel)
  if m.inner.status == 0 || m.inner.status == 1
    return :Optimal
  end
  return :Infeasible
end
getobjval(m::IpoptMathProgModel) = m.inner.obj_val * (m.inner.sense == :Max ? -1 : +1)
getsolution(m::IpoptMathProgModel) = m.inner.x
getconstrsolution(m::IpoptMathProgModel) = m.inner.g
getreducedcosts(m::IpoptMathProgModel) = zeros(m.inner.n)
getconstrduals(m::IpoptMathProgModel) = zeros(m.inner.m)
getrawsolver(m::IpoptMathProgModel) = m.inner


text:
@testset "Aeolian sources" begin
    @testset "$dataset" for dataset in AeolianSources.DATASET_NAMES
        s_A_2D = AeolianSources.load(dataset)
        @testset "$k" for k in keys(s_A_2D)
            v = s_A_2D[k]
            v isa Vector && continue # skip if v is lat/lon
            # Take annual mean if Chien dataset
            v_annual = permutedims(dataset == "Chien" ? dropdims(mean(v, dims=3), dims=3) : v, (2,1))
            # Regrid to OCIM2 grid
            v_regridded = regrid(v_annual, s_A_2D[:lat], s_A_2D[:lon], grd)
            # Paint the top layer
            v_3D = zeros(size(grd)...)
            v_3D[:,:,1] .= ustrip.(upreferred.(v_regridded * u"kg/m^2/s" / grd.δdepth[1]))
            v_vec = v_3D[iwet]
            @test v_vec isa Vector
            @test size(v_vec) == size(iwet)
        end
    end
end

@testset "River sources" begin
    rivers = Rivers.load()
    s = regrid(rivers, grd)
    @test s isa Vector
    @test size(v) == size(iwet)
end

@testset "Groundwater sources" begin
    gws = GroundWaters.load()
    s = regrid(gws, grd)
    @test s isa Vector
    @test size(v) == size(iwet)
end


text:
using Random

n = 256

function apply_seed(seed)
    if seed === missing
        # should result in different numbers across launches
        Random.seed!()
        # XXX: this currently doesn't work, because of the definition in Base,
        #      `seed!(r::MersenneTwister=default_rng())`, which breaks overriding
        #      `default_rng` with a non-MersenneTwister RNG.
    elseif seed !== nothing
        # should result in the same numbers
        Random.seed!(seed)
    elseif seed === nothing
        # should result in different numbers across launches,
        # as determined by the seed set during module loading.
    end
end

@testset "rand($T), seed $seed" for T in (Int32, UInt32, Int64, UInt64, Int128, UInt128,
                                          Float16, Float32, Float64),
                                    seed in (nothing, #=missing,=# 1234)
    # different kernel invocations should get different numbers
    @testset "across launches" begin
        function kernel(A::AbstractArray{T}, seed) where {T}
            apply_seed(seed)
            tid = threadIdx().x
            A[tid] = rand(T)
            return nothing
        end

        a = CUDA.zeros(T, n)
        b = CUDA.zeros(T, n)

        @cuda threads=n kernel(a, seed)
        @cuda threads=n kernel(b, seed)

        if seed === nothing || seed === missing
            @test Array(a) != Array(b)
        else
            @test Array(a) == Array(b)
        end
    end

    # multiple calls to rand should get different numbers
    @testset "across calls" begin
        function kernel(A::AbstractArray{T}, B::AbstractArray{T}, seed) where {T}
            apply_seed(seed)
            tid = threadIdx().x
            A[tid] = rand(T)
            B[tid] = rand(T)
            return nothing
        end

        a = CUDA.zeros(T, n)
        b = CUDA.zeros(T, n)

        @cuda threads=n kernel(a, b, seed)

        @test Array(a) != Array(b)
    end

    # different threads should get different numbers
    @testset "across threads" for active_dim in 1:6
        function kernel(A::AbstractArray{T}, seed) where {T}
            apply_seed(seed)
            id = threadIdx().x*threadIdx().y*threadIdx().z*blockIdx().x*blockIdx().y*blockIdx().z
            A[id] = rand(T)
            return nothing
        end

        tx, ty, tz, bx, by, bz = [dim == active_dim ? 2 : 1 for dim in 1:6]
        a = CUDA.zeros(T, 2)

        @cuda threads=(tx, ty, tz) blocks=(bx, by, bz) kernel(a, seed)

        @test Array(a)[1] != Array(a)[2]
    end
end

@testset "basic randn($T), seed $seed" for T in (Float16, Float32, Float64),
                                           seed in (nothing, #=missing,=# 1234)
    function kernel(A::AbstractArray{T}, seed) where {T}
        apply_seed(seed)
        tid = threadIdx().x
        A[tid] = randn(T)
        return
    end

    a = CUDA.zeros(T, n)
    b = CUDA.zeros(T, n)

    @cuda threads=n kernel(a, seed)
    @cuda threads=n kernel(b, seed)

    if seed === nothing || seed === missing
        @test Array(a) != Array(b)
    else
        @test Array(a) == Array(b)
    end
end

@testset "basic randexp($T), seed $seed" for T in (Float16, Float32, Float64),
                                           seed in (nothing, #=missing,=# 1234)
    function kernel(A::AbstractArray{T}, seed) where {T}
        apply_seed(seed)
        tid = threadIdx().x
        A[tid] = randexp(T)
        return
    end

    a = CUDA.zeros(T, n)
    b = CUDA.zeros(T, n)

    @cuda threads=n kernel(a, seed)
    @cuda threads=n kernel(b, seed)

    if seed === nothing || seed === missing
        @test Array(a) != Array(b)
    else
        @test Array(a) == Array(b)
    end
end


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

using AWS.Compat
using AWS.UUIDs
"""
    BulkPublish()

Initiates a bulk publish of all existing datasets for an Identity Pool to the configured stream. Customers are limited to one successful bulk publish per 24 hours. Bulk publish is an asynchronous request, customers can see the status of the request via the GetBulkPublishDetails operation. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

bulk_publish(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/bulkpublish"; aws_config=aws_config)
bulk_publish(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/bulkpublish", args; aws_config=aws_config)

"""
    DeleteDataset()

Deletes the specific dataset. The dataset will be deleted permanently, and the action can't be undone. Datasets that this dataset was merged with will no longer report the merge. Any subsequent operation on this dataset will result in a ResourceNotFoundException. This API can be called with temporary user credentials provided by Cognito Identity or with developer credentials.

# Required Parameters
- `DatasetName`: A string of up to 128 characters. Allowed characters are a-z, A-Z, 0-9, '_' (underscore), '-' (dash), and '.' (dot).
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

delete_dataset(DatasetName, IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("DELETE", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)"; aws_config=aws_config)
delete_dataset(DatasetName, IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("DELETE", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)", args; aws_config=aws_config)

"""
    DescribeDataset()

Gets meta data about a dataset by identity and dataset name. With Amazon Cognito Sync, each identity has access only to its own data. Thus, the credentials used to make this API call need to have access to the identity data. This API can be called with temporary user credentials provided by Cognito Identity or with developer credentials. You should use Cognito Identity credentials to make this API call.

# Required Parameters
- `DatasetName`: A string of up to 128 characters. Allowed characters are a-z, A-Z, 0-9, '_' (underscore), '-' (dash), and '.' (dot).
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

describe_dataset(DatasetName, IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)"; aws_config=aws_config)
describe_dataset(DatasetName, IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)", args; aws_config=aws_config)

"""
    DescribeIdentityPoolUsage()

Gets usage details (for example, data storage) about a particular identity pool. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

describe_identity_pool_usage(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)"; aws_config=aws_config)
describe_identity_pool_usage(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)", args; aws_config=aws_config)

"""
    DescribeIdentityUsage()

Gets usage information for an identity, including number of datasets and data usage. This API can be called with temporary user credentials provided by Cognito Identity or with developer credentials.

# Required Parameters
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

describe_identity_usage(IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)"; aws_config=aws_config)
describe_identity_usage(IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)", args; aws_config=aws_config)

"""
    GetBulkPublishDetails()

Get the status of the last BulkPublish operation for an identity pool. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

"""

get_bulk_publish_details(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/getBulkPublishDetails"; aws_config=aws_config)
get_bulk_publish_details(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/getBulkPublishDetails", args; aws_config=aws_config)

"""
    GetCognitoEvents()

Gets the events and the corresponding Lambda functions associated with an identity pool. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: The Cognito Identity Pool ID for the request

"""

get_cognito_events(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/events"; aws_config=aws_config)
get_cognito_events(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/events", args; aws_config=aws_config)

"""
    GetIdentityPoolConfiguration()

Gets the configuration settings of an identity pool. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. This is the ID of the pool for which to return a configuration.

"""

get_identity_pool_configuration(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/configuration"; aws_config=aws_config)
get_identity_pool_configuration(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/configuration", args; aws_config=aws_config)

"""
    ListDatasets()

Lists datasets for an identity. With Amazon Cognito Sync, each identity has access only to its own data. Thus, the credentials used to make this API call need to have access to the identity data. ListDatasets can be called with temporary user credentials provided by Cognito Identity or with developer credentials. You should use the Cognito Identity credentials to make this API call.

# Required Parameters
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

# Optional Parameters
- `maxResults`: The maximum number of results to be returned.
- `nextToken`: A pagination token for obtaining the next page of results.
"""

list_datasets(IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets"; aws_config=aws_config)
list_datasets(IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets", args; aws_config=aws_config)

"""
    ListIdentityPoolUsage()

Gets a list of identity pools registered with Cognito. ListIdentityPoolUsage can only be called with developer credentials. You cannot make this API call with the temporary user credentials provided by Cognito Identity.

# Optional Parameters
- `maxResults`: The maximum number of results to be returned.
- `nextToken`: A pagination token for obtaining the next page of results.
"""

list_identity_pool_usage(; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools"; aws_config=aws_config)
list_identity_pool_usage(args::AbstractDict{String, Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools", args; aws_config=aws_config)

"""
    ListRecords()

Gets paginated records, optionally changed after a particular sync count for a dataset and identity. With Amazon Cognito Sync, each identity has access only to its own data. Thus, the credentials used to make this API call need to have access to the identity data. ListRecords can be called with temporary user credentials provided by Cognito Identity or with developer credentials. You should use Cognito Identity credentials to make this API call.

# Required Parameters
- `DatasetName`: A string of up to 128 characters. Allowed characters are a-z, A-Z, 0-9, '_' (underscore), '-' (dash), and '.' (dot).
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.

# Optional Parameters
- `lastSyncCount`: The last server sync count for this record.
- `maxResults`: The maximum number of results to be returned.
- `nextToken`: A pagination token for obtaining the next page of results.
- `syncSessionToken`: A token containing a session ID, identity ID, and expiration.
"""

list_records(DatasetName, IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/records"; aws_config=aws_config)
list_records(DatasetName, IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("GET", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/records", args; aws_config=aws_config)

"""
    RegisterDevice()

Registers a device to receive push sync notifications. This API can only be called with temporary credentials provided by Cognito Identity. You cannot call this API with developer credentials.

# Required Parameters
- `IdentityId`: The unique ID for this identity.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. Here, the ID of the pool that the identity belongs to.
- `Platform`: The SNS platform type (e.g. GCM, SDM, APNS, APNS_SANDBOX).
- `Token`: The push token.

"""

register_device(IdentityId, IdentityPoolId, Platform, Token; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identity/$(IdentityId)/device", Dict{String, Any}("Platform"=>Platform, "Token"=>Token); aws_config=aws_config)
register_device(IdentityId, IdentityPoolId, Platform, Token, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identity/$(IdentityId)/device", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Platform"=>Platform, "Token"=>Token), args)); aws_config=aws_config)

"""
    SetCognitoEvents()

Sets the AWS Lambda function for a given event type for an identity pool. This request only updates the key/value pair specified. Other key/values pairs are not updated. To remove a key value pair, pass a empty value for the particular key. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `Events`: The events to configure
- `IdentityPoolId`: The Cognito Identity Pool to use when configuring Cognito Events

"""

set_cognito_events(Events, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/events", Dict{String, Any}("Events"=>Events); aws_config=aws_config)
set_cognito_events(Events, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/events", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("Events"=>Events), args)); aws_config=aws_config)

"""
    SetIdentityPoolConfiguration()

Sets the necessary configuration for push sync. This API can only be called with developer credentials. You cannot call this API with the temporary user credentials provided by Cognito Identity.

# Required Parameters
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. This is the ID of the pool to modify.

# Optional Parameters
- `CognitoStreams`: Options to apply to this identity pool for Amazon Cognito streams.
- `PushSync`: Options to apply to this identity pool for push synchronization.
"""

set_identity_pool_configuration(IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/configuration"; aws_config=aws_config)
set_identity_pool_configuration(IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/configuration", args; aws_config=aws_config)

"""
    SubscribeToDataset()

Subscribes to receive notifications when a dataset is modified by another device. This API can only be called with temporary credentials provided by Cognito Identity. You cannot call this API with developer credentials.

# Required Parameters
- `DatasetName`: The name of the dataset to subcribe to.
- `DeviceId`: The unique ID generated for this device by Cognito.
- `IdentityId`: Unique ID for this identity.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. The ID of the pool to which the identity belongs.

"""

subscribe_to_dataset(DatasetName, DeviceId, IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/subscriptions/$(DeviceId)"; aws_config=aws_config)
subscribe_to_dataset(DatasetName, DeviceId, IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/subscriptions/$(DeviceId)", args; aws_config=aws_config)

"""
    UnsubscribeFromDataset()

Unsubscribes from receiving notifications when a dataset is modified by another device. This API can only be called with temporary credentials provided by Cognito Identity. You cannot call this API with developer credentials.

# Required Parameters
- `DatasetName`: The name of the dataset from which to unsubcribe.
- `DeviceId`: The unique ID generated for this device by Cognito.
- `IdentityId`: Unique ID for this identity.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. The ID of the pool to which this identity belongs.

"""

unsubscribe_from_dataset(DatasetName, DeviceId, IdentityId, IdentityPoolId; aws_config::AWSConfig=global_aws_config()) = cognito_sync("DELETE", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/subscriptions/$(DeviceId)"; aws_config=aws_config)
unsubscribe_from_dataset(DatasetName, DeviceId, IdentityId, IdentityPoolId, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("DELETE", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)/subscriptions/$(DeviceId)", args; aws_config=aws_config)

"""
    UpdateRecords()

Posts updates to records and adds and deletes records for a dataset and user. The sync count in the record patch is your last known sync count for that record. The server will reject an UpdateRecords request with a ResourceConflictException if you try to patch a record with a new value but a stale sync count. For example, if the sync count on the server is 5 for a key called highScore and you try and submit a new highScore with sync count of 4, the request will be rejected. To obtain the current sync count for a record, call ListRecords. On a successful update of the record, the response returns the new sync count for that record. You should present that sync count the next time you try to update that same record. When the record does not exist, specify the sync count as 0. This API can be called with temporary user credentials provided by Cognito Identity or with developer credentials.

# Required Parameters
- `DatasetName`: A string of up to 128 characters. Allowed characters are a-z, A-Z, 0-9, '_' (underscore), '-' (dash), and '.' (dot).
- `IdentityId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `IdentityPoolId`: A name-spaced GUID (for example, us-east-1:23EC4050-6AEA-7089-A2DD-08002EXAMPLE) created by Amazon Cognito. GUID generation is unique within a region.
- `SyncSessionToken`: The SyncSessionToken returned by a previous call to ListRecords for this dataset and identity.

# Optional Parameters
- `DeviceId`: The unique ID generated for this device by Cognito.
- `RecordPatches`: A list of patch operations.
- `x-amz-Client-Context`: Intended to supply a device ID that will populate the lastModifiedBy field referenced in other methods. The ClientContext field is not yet implemented.
"""

update_records(DatasetName, IdentityId, IdentityPoolId, SyncSessionToken; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)", Dict{String, Any}("SyncSessionToken"=>SyncSessionToken); aws_config=aws_config)
update_records(DatasetName, IdentityId, IdentityPoolId, SyncSessionToken, args::AbstractDict{String, <:Any}; aws_config::AWSConfig=global_aws_config()) = cognito_sync("POST", "/identitypools/$(IdentityPoolId)/identities/$(IdentityId)/datasets/$(DatasetName)", Dict{String, Any}(mergewith(_merge, Dict{String, Any}("SyncSessionToken"=>SyncSessionToken), args)); aws_config=aws_config)


text:
using Documenter, PeriodicSystems
DocMeta.setdocmeta!(PeriodicSystems, :DocTestSetup, :(using PeriodicSystems); recursive=true)

makedocs(
  modules  = [PeriodicSystems],
  sitename = "PeriodicSystems.jl",
  authors  = "Andreas Varga",
  format   = Documenter.HTML(prettyurls = false),
  pages    = [
      "Home"   => "index.md",
      "Library" => [
         "Data Types and Constructors" => [
          "pstypes.md",
          "ps.md"
          ],
          "Basic conversions" =>
          ["psconversions.md",
           "pslifting.md"],
   #      "operations.md",
   #      "order_reduction.md",
          "psanalysis.md",
          "pslyap.md",
   #      "advanced_operations.md",
   #      "model_matching.md"
         ],
     "Utilities" => [
      "pstools.md",
      "psconv.md",
      "slicot.md"
      ],
     "Index" => "makeindex.md"
  ]
)

deploydocs(
  repo = "github.com/andreasvarga/PeriodicSystems.jl.git",
  target = "build",
  devbranch = "master"
)


text:
@testset "Push relabel" begin
    # Construct DiGraph
    flow_graph = DiGraph(8)

    # Load custom dataset
    flow_edges = [
        (1,2,10),(1,3,5),(1,4,15),(2,3,4),(2,5,9),
        (2,6,15),(3,4,4),(3,6,8),(4,7,16),(5,6,15),
        (5,8,10),(6,7,15),(6,8,10),(7,3,6),(7,8,10)
    ]

    capacity_matrix = zeros(Int,8,8)

    for e in flow_edges
        u,v,f = e
        add_edge!(flow_graph,u,v)
        capacity_matrix[u,v] = f
    end
    for g in testdigraphs(flow_graph)
      residual_graph = @inferred(LightGraphs.residual(g))

      # Test enqueue_vertex
      Q = Array{Int,1}()
      excess = [0, 1, 0, 1]
      active = [false, false, true, true]
      @test @inferred(LightGraphs.enqueue_vertex!(Q, 1, active, excess)) == nothing
      @test @inferred(LightGraphs.enqueue_vertex!(Q, 3, active, excess)) == nothing
      @test @inferred(LightGraphs.enqueue_vertex!(Q, 4, active, excess)) == nothing
      @test length(Q) == 0
      @test @inferred(LightGraphs.enqueue_vertex!(Q, 2, active, excess)) == nothing
      @test length(Q) == 1

      # Test push_flow
      Q = Array{Int,1}()
      excess = [15, 1, 1, 0, 0, 0, 0, 0]
      height = [8, 0, 0, 0, 0, 0, 0, 0]
      active = [true, false, false, false, false, false, false, true]
      flow_matrix = zeros(Int, 8, 8)
      @test @inferred(LightGraphs.push_flow!(residual_graph, 1, 2, capacity_matrix, flow_matrix, excess, height, active, Q)) == nothing
      @test length(Q) == 1
      @test flow_matrix[1,2] == 10
      @test @inferred(LightGraphs.push_flow!(residual_graph, 2, 3, capacity_matrix, flow_matrix, excess, height, active, Q)) == nothing
      @test length(Q) == 1
      @test flow_matrix[2,3] == 0

      # Test gap
      Q = Array{Int,1}()
      excess = [15, 1, 1, 0, 0, 0, 0, 0]
      height = [8, 2, 2, 1, 3, 3, 4, 5]
      active = [true, false, false, false, false, false, false, true]
      count  = [0, 1, 2, 2, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0]
      flow_matrix = zeros(Int, 8, 8)

      @test @inferred(LightGraphs.gap!(residual_graph, 1, excess, height, active, count, Q)) == nothing
      @test length(Q) == 2

      # Test relabel
      Q = Array{Int,1}()
      excess = [15, 1, 1, 0, 0, 0, 0, 0]
      height = [8, 1, 1, 1, 1, 1, 1, 0]
      active = [true, false, false, false, false, false, false, true]
      count  = [1, 6, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0]
      flow_matrix = zeros(Int, 8, 8)

      @test @inferred(LightGraphs.relabel!(residual_graph, 2, capacity_matrix, flow_matrix, excess, height, active, count, Q)) == nothing
      @test length(Q) == 1

      # Test discharge
      Q = Array{Int,1}()
      excess = [50, 1, 1, 0, 0, 0, 0, 0]
      height = [8, 0, 0, 0, 0, 0, 0, 0]
      active = [true, false, false, false, false, false, false, true]
      count  = [7, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]
      flow_matrix = zeros(Int, 8, 8)

      @test @inferred(LightGraphs.discharge!(residual_graph, 1, capacity_matrix, flow_matrix, excess, height, active, count, Q)) == nothing
      @test length(Q) == 3

      # Test with default distances
      @test LightGraphs.push_relabel(residual_graph, 1, 8, LightGraphs.DefaultCapacity(residual_graph))[1] == 3

      # Test with capacity matrix
      @test LightGraphs.push_relabel(residual_graph, 1, 8, capacity_matrix)[1] == 28
    end
    # Non regression test added for #448
    M448 =[0 1 0 0 1 1
           1 0 0 0 1 0
           0 0 0 1 0 0
           0 0 0 0 0 0
           1 0 1 0 0 1
           0 0 0 0 1 0]
    g448 = DiGraph(M448)
    @test maximum_flow(g448, 1, 2, M448, algorithm=PushRelabelAlgorithm())[1] == 1
end

max_stars_repo_name: AsileBcd/JohnsonCycles

max_issues_repo_name: AsileBcd/JohnsonCycles

max_forks_repo_name: AsileBcd/JohnsonCycles


text:
using Underscores
using Utils

function parse_file(f)
    @_ map(parse(Int64, _), eachline(f))
end

function load()
    open(joinpath(@__DIR__, "09_encrypted.txt"), "r") do f
        parse_file(f)
    end
end

function search_summands(target, candidates)
    for i in candidates, j in candidates
        if i + j == target
            return i, j
        end
    end
end

function ex1(code, lookback)
    for i in lookback+1:length(code)
        if isnothing(search_summands(code[i], view(code, i-lookback:i-1)))
            return code[i]
        end
    end
end

function ex2(code, nonsum)
    l = length(code)
    for i in eachindex(code)
        s = code[i]
        for j in i+1:l
            s += code[j]
            if s == nonsum
                mi, ma = extrema(view(code, i:j))
                return mi + ma
            elseif s > nonsum
                break
            end
        end
    end
end

code = load()

nonsum = ex1(code, 25)
println("The first number not a sum is: ", nonsum)
println("The code breaking number is ", ex2(code, nonsum))

using Test

@testset "Encoding Error" begin
    example = [
        35,
        20,
        15,
        25,
        47,
        40,
        62,
        55,
        65,
        95,
        102,
        117,
        150,
        182,
        127,
        219,
        299,
        277,
        309,
        576,
    ]

    @testset "example 1" begin
        @test ex1(example, 5) == 127
    end

    @testset "example 2" begin
        @test ex2(example, 127) == 62
    end

    @testset "results" begin
        @test ex1(code, 25) == 31161678
        @test ex2(code, 31161678) == 5453868
    end
end

max_stars_repo_name: CmdQ/AoC


text:
using Documenter
using HITRAN

makedocs(
    modules = [HITRAN],
    authors="Oliver Kliebisch <oliver@kliebisch.net> and contributors",
    repo="https://github.com/tachawkes/hitran.jl/blob/{commit}{path}#L{line}",
    sitename="HITRAN.jl",
    format=Documenter.HTML(;
        prettyurls=get(ENV, "CI", "false") == "true",
        assets=String[],
    ),
    pages=[
        "Introduction & Motivation" => "index.md",
        "Quickstart/Example" => "quickstart.md",
        "Local HITRAN database" => "database.md",
        "Calculating spectra" => "spectra.md"
    ],
)

deploydocs(
    repo = "github.com/TacHawkes/HITRAN.jl.git"
)


text:
"""
    get_int_pos(pos::Real, domain_corner::Real, domain_fac::Real )

Computes the integer position along the PH line.
"""
@inline function get_int_pos(pos::Real, domain_corner::Real, domain_fac::Real )
    return trunc(Int64, ( pos - domain_corner ) * domain_fac)
end

"""
    Lookup tables for peano hilbert keys
"""
const global rottable3 = [36  28  25  27  10  10  25  27;
                          29  11  24  24  37  11  26  26;
                           8   8  25  27  30  38  25  27;
                           9  39  24  24   9  31  26  26;
                          40  24  44  32  40   6  44   6;
                          25   7  33   7  41  41  45  45;
                           4  42   4  46  26  42  34  46;
                          43  43  47  47   5  27   5  35;
                          33  35  36  28  33  35   2   2;
                          32  32  29   3  34  34  37   3;
                          33  35   0   0  33  35  30  38;
                          32  32   1  39  34  34   1  31;
                          24  42  32  46  14  42  14  46;
                          43  43  47  47  25  15  33  15;
                          40  12  44  12  40  26  44  34;
                          13  27  13  35  41  41  45  45;
                          28  41  28  22  38  43  38  22;
                          42  40  23  23  29  39  29  39;
                          41  36  20  36  43  30  20  30;
                          37  31  37  31  42  40  21  21;
                          28  18  28  45  38  18  38  47;
                          19  19  46  44  29  39  29  39;
                          16  36  45  36  16  30  47  30;
                          37  31  37  31  17  17  46  44;
                          12   4   1   3  34  34   1   3;
                           5  35   0   0  13  35   2   2;
                          32  32   1   3   6  14   1   3;
                          33  15   0   0  33   7   2   2;
                          16   0  20   8  16  30  20  30;
                           1  31   9  31  17  17  21  21;
                          28  18  28  22   2  18  10  22;
                          19  19  23  23  29   3  29  11;
                           9  11  12   4   9  11  26  26;
                           8   8   5  27  10  10  13  27;
                           9  11  24  24   9  11   6  14;
                           8   8  25  15  10  10  25   7;
                           0  18   8  22  38  18  38  22;
                          19  19  23  23   1  39   9  39;
                          16  36  20  36  16   2  20  10;
                          37   3  37  11  17  17  21  21;
                           4  17   4  46  14  19  14  46;
                          18  16  47  47   5  15   5  15;
                          17  12  44  12  19  6   44   6;
                          13   7  13   7  18  16  45  45;
                           4  42   4  21  14  42  14  23;
                          43  43  22  20   5  15   5  15;
                          40  12  21  12  40   6  23   6;
                          13   7  13   7  41  41  22  20 ]

const global subpix3 = [0  7  1  6  3  4  2  5 ;
                        7  4  6  5  0  3  1  2 ;
                        4  3  5  2  7  0  6  1 ;
                        3  0  2  1  4  7  5  6 ;
                        1  0  6  7  2  3  5  4 ;
                        0  3  7  4  1  2  6  5 ;
                        3  2  4  5  0  1  7  6 ;
                        2  1  5  6  3  0  4  7 ;
                        6  1  7  0  5  2  4  3 ;
                        1  2  0  3  6  5  7  4 ;
                        2  5  3  4  1  6  0  7 ;
                        5  6  4  7  2  1  3  0 ;
                        7  6  0  1  4  5  3  2 ;
                        6  5  1  2  7  4  0  3 ;
                        5  4  2  3  6  7  1  0 ;
                        4  7  3  0  5  6  2  1 ;
                        6  7  5  4  1  0  2  3 ;
                        7  0  4  3  6  1  5  2 ;
                        0  1  3  2  7  6  4  5 ;
                        1  6  2  5  0  7  3  4 ;
                        2  3  1  0  5  4  6  7 ;
                        3  4  0  7  2  5  1  6 ;
                        4  5  7  6  3  2  0  1 ;
                        5  2  6  1  4  3  7  0 ;
                        7  0  6  1  4  3  5  2 ;
                        0  3  1  2  7  4  6  5 ;
                        3  4  2  5  0  7  1  6 ;
                        4  7  5  6  3  0  2  1 ;
                        6  7  1  0  5  4  2  3 ;
                        7  4  0  3  6  5  1  2 ;
                        4  5  3  2  7  6  0  1 ;
                        5  6  2  1  4  7  3  0 ;
                        1  6  0  7  2  5  3  4 ;
                        6  5  7  4  1  2  0  3 ;
                        5  2  4  3  6  1  7  0 ;
                        2  1  3  0  5  6  4  7 ;
                        0  1  7  6  3  2  4  5 ;
                        1  2  6  5  0  3  7  4 ;
                        2  3  5  4  1  0  6  7 ;
                        3  0  4  7  2  1  5  6 ;
                        1  0  2  3  6  7  5  4 ;
                        0  7  3  4  1  6  2  5 ;
                        7  6  4  5  0  1  3  2 ;
                        6  1  5  2  7  0  4  3 ;
                        5  4  6  7  2  3  1  0 ;
                        4  3  7  0  5  2  6  1 ;
                        3  2  0  1  4  5  7  6 ;
                        2  5  1  6  3  4  0  7 ]

"""
    peano_hilbert_key(bits::Integer, x::Integer, y::Integer, z::Integer)

Computes a Peano-Hilbert key for an integer triplet (x,y,z) with x,y,z typically in the range
between 0 and 2^bits-1. Values outside this range can occur when reading across the
periodic box borders.
"""
function peano_hilbert_key(bits::Integer, x::Integer, y::Integer, z::Integer)

    rotation = 1
    key = 0
    mask = 1 << (bits-1)

    @inbounds for _ ∈ 1:bits

        pix = (
                (( x & mask) > 0 ? 4 : 0)  +
                (( y & mask) > 0 ? 2 : 0)  +
                (( z & mask) > 0 ? 1 : 0)
              ) + 1

        key = key << 3

        key = key | subpix3[rotation, pix]

        rotation = rottable3[rotation, pix] + 1

        mask = mask >> 1

    end

    return key
end

"""
    get_keylist(h_key::KeyHeader, x0::Array{<:Real}, x1::Array{<:Real})

Get all Peano-Hilbert keys for domain defined by the corner points `x0` and `x1`.
"""
function get_keylist(h_key::KeyHeader, x0::Array{T}, x1::Array{T}) where T

    ix0 = Vector{Int}(undef, 3)
    ix1 = Vector{Int}(undef, 3)
    dix = Vector{Int}(undef, 3)

    nkeys = 1

    @inbounds for i = 1:3
        ix0[i] = get_int_pos( x0[i], h_key.domain_corners[i], h_key.domain_fac )
        ix1[i] = get_int_pos( x1[i], h_key.domain_corners[i], h_key.domain_fac )
        dix[i] = ix1[i] - ix0[i] + 1
        nkeys *= dix[i]
    end

    keylist = Vector{Int}(undef, nkeys)

    i = 1
    @inbounds for ix = ix0[1]:ix1[1], iy = ix0[2]:ix1[2], iz = ix0[3]:ix1[3]
        keylist[i] = peano_hilbert_key(h_key.bits, ix, iy, iz)
        i += 1
    end

    return sort!(keylist)
end


text:
####################################
## check for dataarray with idx ##
####################################

function chkIdx(idx)
    ## check for Array and entry type
    ## Array is already checked through Timedata definition itself
    if !isa(idx, Array)
        throw(TypeError(:chkIdx,
                        "passed time index: must be given as array\n",
                        Array,
                        typeof(idx)))
    end
    if !issubtype(eltype(idx), Union(Integer, Float64, Date, DateTime))
        throw(TypeError(:chkIdx,
                        "passed time index: entries must be either
Integer, Float64 or TimeType\n",
                        Union(Integer, Date, DateTime),
                        eltype(idx)))
    end
end

###########################################
## check for numeric values in dataframe ##
###########################################

function isNumericType(typ::Type)
    # numeric and NAtype allowed, boolean not
    valid = (issubtype(typ, Union(Number, NAtype))) &
               (!issubtype(typ, Bool))
    return valid
end

function chkNumDf(df::DataFrame)
    ## check for numeric values or NAs

    n = ncol(df)
    errMsg = "all columns must be numeric for conversion"
    typs = eltypes(df)
    for ii=1:n
        ## check for numeric values
        if !isNumericType(typs[ii])
            throw(ArgumentError(errMsg))
        end
    end
end

############################################
## check for numeric values only - no NAs ##
############################################

function chkNum(df::DataFrame)
    ## check for numeric values, no NAs allowed

    for ii=1:size(df, 2)
        if any(isna(df[ii]))
            throw(ArgumentError("no NAs allowed in Timematr"))
        end
    end
    chkNumDf(df)
end

###########################################
## check for values within unit interval ##
###########################################

function chkUnit(df::DataFrame)
    vals = convert(Array, df)
    if (any(vals .< 0) | any(vals .> 1))
        throw(ArgumentError("values must be inside of unit interval"))
    end
end


text:
@testset "Measurement" begin

    @testset "Multiple Emitters" begin
        manifold = IdealManifold()
        system = GPSL1()

        sat1 = ConstantDopplerSatellite(
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sat2 = ConstantDopplerSatellite(
            system,
            2,
            carrier_doppler = 500.0Hz,
            carrier_phase = π / 2,
            code_phase = 50.0,
            cn0 = 45dBHz
        )
        sats = (sat1,sat2)
        receiver = @inferred Receiver(2.5e6Hz)

        signal = StructArray{ComplexF32}(undef, 2500)
        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats = @inferred get_measurement!(
            signal,
            receiver,
            sats,
            manifold,
            rng
        )

        carrier = cis.(2π .* 1000.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 1000 / 1540) ./ 2.5e6 .+ 100, 1) .*
            sqrt(10^(45 / 10))
        carrier = cis.(2π .* 500.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.25)
        reference_signal += carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 500 / 1540) ./ 2.5e6 .+ 50, 2) .*
            sqrt(10^(45 / 10))
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF64, 2500) .* sqrt(2.5e6)
        @test measurement ≈ reference_signal

        @test next_receiver == GNSSSimulator.propagate(receiver, 2500, rng)
        @test next_sats == GNSSSimulator.propagate.(sats, 2500, 0.0Hz, 2.5e6Hz, 1.0/Hz, Ref(rng))

        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats = @inferred get_measurement(
            2500,
            receiver,
            sats,
            manifold,
            rng
        )
        @test measurement ≈ reference_signal

        sat1f64 = ConstantDopplerSatellite(
            Float64,
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sat2f64 = ConstantDopplerSatellite(
            Float64,
            system,
            2,
            carrier_doppler = 500.0Hz,
            carrier_phase = π / 2,
            code_phase = 50.0,
            cn0 = 45dBHz
        )
        satsf64 = (sat1f64,sat2f64)
        receiver = @inferred Receiver(2.5e6Hz)

        signalf64 = StructArray{ComplexF64}(undef, 2500)
        rng = MersenneTwister(1234)
        measurementf64, next_receiverf64, next_satsf64 = @inferred get_measurement!(
            signalf64,
            receiver,
            satsf64,
            manifold,
            rng
        )
        @test measurementf64 ≈ reference_signal

        rng = MersenneTwister(1234)
        measurementf64, next_receiver, next_sats = @inferred get_measurement(
            2500,
            receiver,
            satsf64,
            manifold,
            rng
        )
        @test measurementf64 ≈ reference_signal
    end

    @testset "Non existing Emitter" begin
        manifold = IdealManifold()
        system = GPSL1()

        sat1 = ConstantDopplerSatellite(
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sat2 = ConstantDopplerSatellite(
            system,
            2,
            carrier_doppler = 500.0Hz,
            carrier_phase = π / 2,
            code_phase = 50.0,
            cn0 = 45dBHz,
            exists = false
        )
        sats = (sat1,sat2)
        receiver = @inferred Receiver(2.5e6Hz)

        signal = StructArray{ComplexF32}(undef, 2500)
        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats = @inferred get_measurement!(
            signal,
            receiver,
            sats,
            manifold,
            rng
        )

        carrier = cis.(2π .* 1000.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 1000 / 1540) ./ 2.5e6 .+ 100, 1) .*
            sqrt(10^(45 / 10))
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF32, 2500) .* sqrt(2.5e6)
        @test measurement ≈ reference_signal
    end

    @testset "Signal continuity" begin
        manifold = IdealManifold()
        system = GPSL1()

        sat1 = ConstantDopplerSatellite(
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sats = (sat1,)
        receiver = @inferred Receiver(2.5e6Hz)

        rng = MersenneTwister(1234)
        measurement1, next_receiver, next_sats = @inferred get_measurement(
            2500,
            receiver,
            sats,
            manifold,
            rng
        )

        @test next_receiver == GNSSSimulator.propagate(receiver, 2500, rng)
        @test next_sats == GNSSSimulator.propagate.(sats, 2500, 0.0Hz, 2.5e6Hz, 1.0/Hz, Ref(rng))

        measurement2, next_receiver, next_sats = @inferred get_measurement(
            2500,
            next_receiver,
            next_sats,
            manifold,
            rng
        )
        measurement = [measurement1; measurement2]

        carrier = cis.(2π .* 1000.0Hz ./ 2.5e6Hz .* (0:4999) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(system, (0:4999) .* (1023e3 + 1000 / 1540) ./ 2.5e6 .+ 100, 1) .*
            sqrt(10^(45 / 10))
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF32, 5000) .* sqrt(2.5e6)

        @test measurement.im ≈ imag.(reference_signal)
        # Slightly not contineous because of Float and Fixed point representation?
        @test measurement.re[1:2500] ≈ real.(reference_signal[1:2500])
        @test measurement.re[2502:4917] ≈ real.(reference_signal[2502:4917])
        @test measurement.re[4919:end] ≈ real.(reference_signal[4919:end])
    end

    @testset "Multiple antennas" begin
        manifold = IdealManifold(
            1575420e3,
            0.1904 / 4 * SVector(SVector(1, 1, 0), SVector(-1, 1, 0)),
        )
        system = GPSL1()

        sat1 = ConstantDopplerSatellite(
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sats = (sat1,)
        gpmc = @SMatrix randn(ComplexF64, 2, 2)
        receiver = @inferred Receiver(
            2.5e6Hz,
            gain_phase_mism_crosstalk = gpmc
        )

        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats = @inferred get_measurement(
            2500,
            receiver,
            sats,
            manifold,
            rng
        )
        carrier = cis.(2π .* 1000.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 1000 / 1540) ./ 2.5e6 .+ 100, 1) .*
            sqrt(10^(45 / 10)) .* transpose([1.0, 1.0])
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF32, 2500, 2) .* sqrt(2.5e6)
        reference_signal = reference_signal * transpose(gpmc)
        @test measurement ≈ reference_signal
    end

    @testset "Two emitter types" begin
        manifold = IdealManifold()
        gpsl1 = GPSL1()
        gpsl5 = GPSL5()

        sat1 = ConstantDopplerSatellite(
            gpsl1,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sat2 = ConstantDopplerSatellite(
            gpsl5,
            1,
            carrier_doppler = 500.0Hz,
            carrier_phase = π / 2,
            code_phase = 50.0,
            cn0 = 45dBHz
        )
        sats = (sat1,sat2)
        receiver = @inferred Receiver(25e6Hz)

        jammer = CWJammer(1, 10dB, doppler = 100.0Hz, phase = π / 4)
        jammers = (jammer, )

        num_samples = 25000

        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats, next_jammers = @inferred get_measurement(
            num_samples,
            receiver,
            sats,
            jammers,
            manifold,
            rng
        )

        carrier = cis.(2π .* 1000.0Hz ./ 25e6Hz .* (0:num_samples-1) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(gpsl1, (0:num_samples-1) .* (1023e3 + 1000 / 1540) ./ 25e6 .+ 100, 1) .*
            sqrt(10^(45 / 10))
        carrier = cis.(2π .* 500.0Hz ./ 25e6Hz .* (0:num_samples-1) .+ 2π .* 0.25)
        reference_signal += carrier .*
            get_code.(gpsl5, (0:num_samples-1) .* (1023e4 + 500 / 115) ./ 25e6 .+ 50, 1) .*
            sqrt(10^(45 / 10))
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF64, num_samples) .* sqrt(25e6)
        carrier = cis.(2π .* 100.0Hz ./ 25e6Hz .* (0:num_samples-1) .+ 2π .* 0.125)
        reference_signal += carrier .* sqrt(10^(10 / 10) * 25e6)
        @test measurement ≈ reference_signal

    end

    @testset "Three emitter types" begin
        manifold = IdealManifold()
        system = GPSL1()

        sat1 = ConstantDopplerSatellite(
            system,
            1,
            carrier_doppler = 1000.0Hz,
            carrier_phase = π / 2,
            code_phase = 100.0,
            cn0 = 45dBHz
        )
        sats = (sat1, )
        receiver = @inferred Receiver(2.5e6Hz)

        jammer = CWJammer(1, 10dB, doppler = 100.0Hz, phase = π / 4)
        jammers = (jammer, )

        si = ConstantDopplerStructuralInterference(
            sat1,
            -3dB,
            added_carrier_doppler = 10Hz,
            added_carrier_phase = π / 8,
            added_code_phase = 10.0
        )

        sis = (si, )

        rng = MersenneTwister(1234)
        measurement, next_receiver, next_sats, next_jammers, next_sis = @inferred get_measurement(
            2500,
            receiver,
            sats,
            jammers,
            sis,
            manifold,
            rng
        )

        carrier = cis.(2π .* 1000.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.25)
        reference_signal = carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 1000 / 1540) ./ 2.5e6 .+ 100, 1) .*
            sqrt(10^(45 / 10))
        carrier = cis.(2π .* 1010.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.3125)
        reference_signal += carrier .*
            get_code.(system, (0:2499) .* (1023e3 + 1010 / 1540) ./ 2.5e6 .+ 110, 1) .*
            sqrt(10^(42 / 10))
        rng = MersenneTwister(1234)
        reference_signal += randn(rng, ComplexF64, 2500) .* sqrt(2.5e6)
        carrier = cis.(2π .* 100.0Hz ./ 2.5e6Hz .* (0:2499) .+ 2π .* 0.125)
        reference_signal += carrier .* sqrt(10^(10 / 10) * 2.5e6)
        @test measurement ≈ reference_signal
    end
end


text:
"""
Construct model for HydroGen with FixedOutput Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, FixedOutput},
    ::Type{S},
) where {H <: PSY.HydroGen, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    nodal_expression!(psi_container, devices, S)

    return
end

"""
Construct model for HydroGen with RunOfRiver Dispatch Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {
    H <: PSY.HydroGen,
    D <: AbstractHydroDispatchFormulation,
    S <: PM.AbstractPowerModel,
}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    add_constraints!(
        psi_container,
        RangeConstraint,
        ReactivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with RunOfRiver Dispatch Formulation
with only Active Power.
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {
    H <: PSY.HydroGen,
    D <: AbstractHydroDispatchFormulation,
    S <: PM.AbstractActivePowerModel,
}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirBudget Dispatch Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchReservoirBudget},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)

    #Energy Budget Constraint
    energy_budget_constraints!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirBudget Dispatch Formulation
with only Active Power.
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchReservoirBudget},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractActivePowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)

    #Energy Budget Constraint
    energy_budget_constraints!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirStorage Dispatch Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchReservoirStorage},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)
    add_variables!(psi_container, EnergyVariable, devices)
    add_variables!(psi_container, SpillageVariable, devices)

    #Initial Conditions
    storage_energy_init(psi_container, devices)
    #Energy Balance Constraint
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirStorage Dispatch Formulation
with only Active Power
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchReservoirStorage},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractActivePowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, EnergyVariable, devices)
    add_variables!(psi_container, SpillageVariable, devices)

    #Initial Conditions
    storage_energy_init(psi_container, devices)
    #Energy Balance Constraint
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with RunOfRiver Commitment Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {H <: PSY.HydroGen, D <: HydroCommitmentRunOfRiver, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ReactivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    commit_hydro_active_power_ub!(psi_container, devices, model, get_feedforward(model))
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with RunOfRiver Commitment Formulation
with only Active Power.
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {
    H <: PSY.HydroGen,
    D <: HydroCommitmentRunOfRiver,
    S <: PM.AbstractActivePowerModel,
}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    commit_hydro_active_power_ub!(psi_container, devices, model, get_feedforward(model))
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirBudget Commitment Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {H <: PSY.HydroGen, D <: HydroCommitmentReservoirBudget, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    add_constraints!(
        psi_container,
        RangeConstraint,
        ReactivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    #Energy Budget Constraint
    energy_budget_constraints!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirBudget Commitment Formulation
with only Active Power.
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, D},
    ::Type{S},
) where {
    H <: PSY.HydroGen,
    D <: HydroCommitmentReservoirBudget,
    S <: PM.AbstractActivePowerModel,
}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    #Energy Budget Constraint
    energy_budget_constraints!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirStorage Commitment Formulation
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroCommitmentReservoirStorage},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractPowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, ReactivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)
    add_variables!(psi_container, EnergyVariable, devices)
    add_variables!(psi_container, SpillageVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    add_constraints!(
        psi_container,
        RangeConstraint,
        ReactivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )

    #Initial Conditions
    storage_energy_init(psi_container, devices)
    #Energy Balance Constraint
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroGen with ReservoirStorage Dispatch Formulation
with only Active Power
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroCommitmentReservoirStorage},
    ::Type{S},
) where {H <: PSY.HydroEnergyReservoir, S <: PM.AbstractActivePowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerVariable, devices)
    add_variables!(psi_container, OnVariable, devices)
    add_variables!(psi_container, EnergyVariable, devices)
    add_variables!(psi_container, SpillageVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )

    #Initial Conditions
    storage_energy_init(psi_container, devices)
    #Energy Balance Constraint
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))
    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function!(psi_container, devices, model, S, nothing)

    return
end

"""
Construct model for HydroPumpedStorage with PumpedStorage Dispatch Formulation
with only Active Power
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchPumpedStorage},
    ::Type{S},
) where {H <: PSY.HydroPumpedStorage, S <: PM.AbstractActivePowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerInVariable, devices)
    add_variables!(psi_container, ActivePowerOutVariable, devices)
    add_variables!(psi_container, EnergyVariableUp, devices)
    add_variables!(psi_container, EnergyVariableDown, devices)
    add_variables!(psi_container, SpillageVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerOutVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerInVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )

    #Initial Conditions
    storage_energy_init(psi_container, devices)

    # Energy Balanace limits
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function(psi_container, devices, HydroDispatchReservoirBudget, S)

    return
end

"""
Construct model for HydroPumpedStorage with PumpedStorage Dispatch Formulation with
reservation constraint with only Active Power
"""
function construct_device!(
    psi_container::PSIContainer,
    sys::PSY.System,
    model::DeviceModel{H, HydroDispatchPumpedStoragewReservation},
    ::Type{S},
) where {H <: PSY.HydroPumpedStorage, S <: PM.AbstractActivePowerModel}
    devices = get_available_components(H, sys)

    if !validate_available_devices(H, devices)
        return
    end

    #Variables
    add_variables!(psi_container, ActivePowerInVariable, devices)
    add_variables!(psi_container, ActivePowerOutVariable, devices)
    add_variables!(psi_container, EnergyVariableUp, devices)
    add_variables!(psi_container, EnergyVariableDown, devices)
    add_variables!(psi_container, SpillageVariable, devices)
    add_variables!(psi_container, ReserveVariable, devices)

    #Constraints
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerOutVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )
    add_constraints!(
        psi_container,
        RangeConstraint,
        ActivePowerInVariable,
        devices,
        model,
        S,
        get_feedforward(model),
    )

    #Initial Conditions
    storage_energy_init(psi_container, devices)

    # Energy Balanace limits
    energy_balance_constraint!(psi_container, devices, model, S, get_feedforward(model))

    feedforward!(psi_container, devices, model, get_feedforward(model))

    #Cost Function
    cost_function(psi_container, devices, HydroDispatchReservoirBudget, S)

    return
end


text:
module Indexing

import Base: getindex, setindex!, view, IndexStyle, @propagate_inbounds, keys, haskey,
             parent, length, iterate, axes, size

export getindices, setindices!, ViewArray, ViewVector, ViewMatrix, ViewDict

_valtype(a::AbstractDict) = Base.valtype(a)
_valtype(a) = Base.eltype(a)

include("getindices.jl")
include("setindices.jl")
include("view.jl")

# TODO:
#
# * Deal better with bounds checking and @inbounds
# * Accelerate for Base containers where necessary?

end # module


text:
"""
    function assign_lines_index(lines)

Finds the binding lines of a power system which has been solved by OPF, and reorders them in
accordance with their assigned index.

#Arguments
- `lines`:                      lines as a direct output from
                                lines = get_components(Line, sys), see PowerSystems.jl

# Throws
- `ERROR`:                      lines has no field name. The argument must be a direct
                                output from lines = get_components(Line, sys),
                                see PowerSystems.jl

"""
function assign_lines_index(lines)
    # Get lines names and sort them
    all_lines = Array{Any}(undef,(length(lines),2))
    line_count = 0
    for line in lines
        line_count = line_count + 1
        all_lines[line_count,1] = line.name
    end

    all_lines_sorted = Array{Any}(undef,(length(lines), 2))
    all_lines_sorted[:, 1] = sort(all_lines[:, 1])
    all_lines_sorted[:, 2] = [1:length(lines);]
    # Assign Index number to unsorted lines
    for ul = 1 : length(lines) #unsorted lines
        for sl = 1 : length(lines) #sorted lines
            if all_lines_sorted[sl, 1] == all_lines[ul, 1]
                all_lines[ul, 2] = all_lines_sorted[sl, 2]
            end
        end
    end

    return all_lines, all_lines_sorted
end

max_stars_repo_name: invenia/MaximumGeneratorProfitTCRDD

max_issues_repo_name: invenia/MaximumGeneratorProfitTCRDD

max_forks_repo_name: invenia/MaximumGeneratorProfitTCRDD


text:
@reexport module DayAhead

# Standard library
using Statistics
using Printf
using Dates

using RecipesBase
using Distributions
using Parameters
using StochasticPrograms
using Flux
using HydroModels
using HydroModels: AbstractHydroModel, StochasticHydroModel, River, Plant, Area, Scenario, Forecaster, forecast, Q̄, μ, %

import Base: show
import HydroModels: modelindices

using Plots: font, text, Shape, RGB

export
    DayAheadData,
    NordPoolDayAheadData,
    DayAheadScenario,
    DayAheadScenarios,
    RecurrentDayAheadSampler,
    DayAheadModel,
    strategy,
    singleorder,
    singleorders,
    blockorders,
    independent,
    dependent

include("data.jl")
include("scenarios.jl")
include("indices.jl")
include("model.jl")
include("orderstrategy.jl")

end


text:
using HybridSystems
using Polyhedra
using SwitchOnSafety
using Test

include("solvers.jl")

include("jsr.jl")
include("invariant.jl")

include("../examples/run_examples.jl")


text:
export load_images, plot_image

using FITSIO
using AstroLib
using Plots

@doc raw"""
    load_images(path::String)

Load images from the `path` directory and its time-stamp. There are two different types
of time-stamp loading, the first one, if "JD" is header entry, so it is taken by that, else (for example, the well-known software CCDSoft, provided by SBIG company) it converts a date-obs into JD. If you are looking for light curve TOM, you have to apply the heliocentric correction.

The heliocentric correction of the date-obs has not been implemented yet.
"""
function load_images(path::String)

    file_list = readdir(path)
    img_stack = Array{Array{UInt16, 2}, 1}(undef, 0)
    times = Array{Float64, 1}(undef, 0)

    for i in file_list
        fit = FITS(path*i)
        img = read(fit[1])

        try
            time_obs = read_header(fit[1])["JD"]
            times = append!(times, time_obs)
        catch e
            jd = read_header(fit[1])["DATE-OBS"]
            time_obs = jdcnv(jd)
            times = append!(times, time_obs)
        end

        img_stack = append!(img_stack, [img])

    end

    return img_stack, times

end

@doc raw"""
    plot_image(img::Array{UInt16,2}, color_lims::Tuple)

A function to easily plot images. It is based on Plots.
"""
function plot_image(img::Array{UInt16,2}, color_lims::Tuple)
    gradient=ColorGradient([:black, :white])
    heatmap(transpose(img), c=gradient, clims=color_lims)
end


text:
println("Starting precompilation")
using OrdinaryDiffEq
using Flux, DiffEqFlux, GalacticOptim
using Test
using Distributions
using Distances
using Zygote
using DistributionsAD
using LinearAlgebra
using Random
println("Starting tests")

Random.seed!(1999)

## callback to be used by all tests
function cb(p, l)
    @show l
    false
end

###
# test for default base distribution and deterministic trace CNF
###

nn = Chain(Dense(1, 1, tanh))
data_dist = Beta(7, 7)
data_train = Float32.(rand(data_dist, 1, 100))
tspan = (0.0f0, 1.0f0)
cnf_test = FFJORD(nn, tspan, Tsit5())

function loss_adjoint(θ)
    logpx = cnf_test(data_train, θ; monte_carlo=false)[1]
    loss = -mean(logpx)
end

optfunc = GalacticOptim.OptimizationFunction((x, p) -> loss_adjoint(x), GalacticOptim.AutoZygote())
optprob = GalacticOptim.OptimizationProblem(optfunc, cnf_test.p)
res = GalacticOptim.solve(optprob, ADAM(0.1), cb=cb, maxiters=100)

θopt = res.minimizer
data_validate = Float32.(rand(data_dist, 1, 100))
actual_pdf = pdf.(data_dist, data_validate)
# use direct trace calculation for predictions
learned_pdf = exp.(cnf_test(data_validate, θopt; monte_carlo=false)[1])

@test totalvariation(learned_pdf, actual_pdf) / size(data_validate, 2) < 0.35

###
# test for alternative base distribution and deterministic trace CNF
###

nn = Chain(Dense(1, 3, tanh), Dense(3, 1, tanh))
data_dist = Normal(6.0, 0.7)
data_train = Float32.(rand(data_dist, 1, 100))
tspan = (0.0f0, 1.0f0)
cnf_test = FFJORD(nn, tspan, Tsit5(); basedist=MvNormal([0.0f0], Diagonal([4.0f0])))

function loss_adjoint(θ)
    logpx = cnf_test(data_train, θ; monte_carlo=false)[1]
    loss = -mean(logpx)
end

optfunc = GalacticOptim.OptimizationFunction((x, p) -> loss_adjoint(x), GalacticOptim.AutoZygote())
optprob = GalacticOptim.OptimizationProblem(optfunc, 0.01f0 .* cnf_test.p)
res = GalacticOptim.solve(optprob, ADAM(0.1), cb=cb, maxiters=100)

θopt = res.minimizer
data_validate = Float32.(rand(data_dist, 1, 100))
actual_pdf = pdf.(data_dist, data_validate)
learned_pdf = exp.(cnf_test(data_validate, θopt; monte_carlo=false)[1])

@test totalvariation(learned_pdf, actual_pdf) / size(data_validate, 2) < 0.25

###
# test for multivariate distribution and deterministic trace CNF
###

nn = Chain(Dense(2, 2, tanh))
μ = ones(2)
Σ = Diagonal([7.0, 7.0])
data_dist = MvNormal(μ, Σ)
data_train = Float32.(rand(data_dist, 100))
tspan = (0.0f0, 1.0f0)
cnf_test = FFJORD(nn, tspan, Tsit5())

function loss_adjoint(θ)
    logpx = cnf_test(data_train, θ; monte_carlo=false)[1]
    loss = -mean(logpx)
end

optfunc = GalacticOptim.OptimizationFunction((x, p) -> loss_adjoint(x), GalacticOptim.AutoZygote())
optprob = GalacticOptim.OptimizationProblem(optfunc, 0.01f0 .* cnf_test.p)
res = GalacticOptim.solve(optprob, ADAM(0.1), cb=cb, maxiters=300)

θopt = res.minimizer
data_validate = Float32.(rand(data_dist, 100))
actual_pdf = pdf(data_dist, data_validate)
learned_pdf = exp.(cnf_test(data_validate, θopt; monte_carlo=false)[1])

@test totalvariation(learned_pdf, actual_pdf) / size(data_validate, 2) < 0.25

###
# test for default multivariate distribution and FFJORD with regularizers
###

nn = Chain(Dense(1, 1, tanh))
data_dist = Beta(7, 7)
data_train = Float32.(rand(data_dist, 1, 100))
tspan = (0.0f0, 1.0f0)
ffjord_test = FFJORD(nn, tspan, Tsit5())

function loss_adjoint(θ)
    logpx, λ₁, λ₂ = ffjord_test(data_train, θ, true)
    return mean(@. -logpx + 0.1 * λ₁ + λ₂)
end

optfunc = GalacticOptim.OptimizationFunction((x, p) -> loss_adjoint(x), GalacticOptim.AutoZygote())
optprob = GalacticOptim.OptimizationProblem(optfunc, 0.01f0 .* ffjord_test.p)
res = GalacticOptim.solve(optprob, ADAM(0.1), cb=cb, maxiters=300)

θopt = res.minimizer
data_validate = Float32.(rand(data_dist, 1, 100))
actual_pdf = pdf.(data_dist, data_validate)
learned_pdf = exp.(ffjord_test(data_validate, θopt; monte_carlo=false)[1])

@test totalvariation(learned_pdf, actual_pdf) / size(data_validate, 2) < 0.40


text:
using Test

@testset "Core" begin
  include("core/Core.jl")
end

@testset "Theories" begin
  include("theories/Theories.jl")
end

@testset "CategoricalAlgebra" begin
  include("categorical_algebra/CategoricalAlgebra.jl")
end

@testset "Graphs" begin
  include("graphs/Graphs.jl")
end

@testset "WiringDiagrams" begin
  include("wiring_diagrams/WiringDiagrams.jl")
end

@testset "Graphics" begin
  include("graphics/Graphics.jl")
end

@testset "Programs" begin
  include("programs/Programs.jl")
end

@testset "LinearAlgebra" begin
  include("linear_algebra/LinearAlgebra.jl")
end


text:
module LSMH

using Combinatorics, BioSequences

export
    murmur32,
    xxhash32,
    minhash, signature, lsh,
    fasta_records, fasta_signatures,
    AbstractSignature, Signature, identifier, hashband,
    HashTable, band, hashtable,
    candidates, jaccard

include("hashing.jl")
include("murmurhash3.jl")
include("xxhash.jl")
include("proteins.jl")

end  # module


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

using FEMBase, MortarContact2DAD, Test
using MortarContact2DAD: get_slave_dofs, get_master_dofs
using MortarContact2DAD: project_from_master_to_slave_ad, project_from_slave_to_master_ad

X = Dict(
    1 => [0.0, 0.0],
    2 => [1.0, 0.0],
    3 => [0.0, 1.0],
    4 => [1.0, 1.0])

u = Dict(
    1 => [0.0, 0.0],
    2 => [0.0, 0.0],
    3 => [0.0, 0.0],
    4 => [0.0, 0.0])

slave = Element(Seg2, [1, 2])
master = Element(Seg2, [3, 4])
update!([slave, master], "geometry", X)
update!([slave, master], "displacement", u)
update!(slave, "normal", ([0.0, 1.0], [0.0, 1.0]))
problem = Problem(Mortar2DAD, "test problem", 2, "displacement")
add_slave_elements!(problem, [slave])
add_master_elements!(problem, [master])
problem.assembly.u = zeros(8)
problem.assembly.la = zeros(8)
assemble!(problem, 0.0)
C1_expected = 1/6*[
               2.0 0.0 1.0 0.0 -2.0 0.0 -1.0 0.0
               0.0 2.0 0.0 1.0 0.0 -2.0 0.0 -1.0
               1.0 0.0 2.0 0.0 -1.0 0.0 -2.0 0.0
               0.0 1.0 0.0 2.0 0.0 -1.0 0.0 -2.0]
C1 = problem.assembly.C1
C2 = problem.assembly.C2
@test isapprox(C1, C2)
@test isapprox(C1, C1_expected)

@test_throws Exception add_elements!(problem, [slave])
@test get_slave_dofs(problem) == [1, 2, 3, 4]
@test get_master_dofs(problem) == [5, 6, 7, 8]
x1 = slave("geometry", 0.0)
n1 = slave("normal", 0.0)
xm = 1/3*(X[3]+X[4])
xs = 1/3*(X[1]+X[2])
@test_throws ErrorException project_from_master_to_slave_ad(slave, x1, n1, xm, 0.0; tol=0.0, max_iterations=1)
x1 = [0.5, 0.0]
n1 = [0.0, 1.0]
x2 = master("geometry", 0.0)
@test_throws ErrorException project_from_slave_to_master_ad(master, x1, n1, x2, 0.0; tol=0.0, max_iterations=1)


text:
using LinearAlgebraicRepresentation, ViewerGL, SparseArrays
Lar = LinearAlgebraicRepresentation; GL = ViewerGL
using QHull

V,CV = Lar.rod()([4,1])
ch = QHull.chull(convert(Lar.Points,V'))
FV = ch.simplices
EV = Lar.simplexFacets(FV)
GL.VIEW([ GL.GLGrid(V,EV) ]);

cyl = Lar.Struct([ Lar.t(0,0,-1.5), (V,FV,EV) ])
tris = Lar.Struct([ cyl,
                    Lar.Struct([Lar.r(pi/2,0,0), cyl]),
                    Lar.Struct([Lar.r(0,pi/2,0), cyl])
                    ])
V,FV,EV = Lar.struct2lar(tris)
GL.VIEW([ GL.GLGrid(V,FV) ]);

VV = [[k] for k=1:size(V,2)]
GL.VIEW( GL.numbering(.2)((V,[VV,EV,FV]), GL.COLORS[4]) );

cop_EV = convert(Lar.ChainOp, Lar.coboundary_0(EV::Lar.Cells));
cop_FE = Lar.coboundary_1(V, FV::Lar.Cells, EV::Lar.Cells);
W = convert(Lar.Points, V');
V, copEV, copFE, copCF = Lar.Arrangement.spatial_arrangement( W, cop_EV, cop_FE)

W = convert(Lar.Points, V')
V,CVs,FVs,EVs = Lar.pols2tria(W, copEV, copFE, copCF)

GL.VIEW(GL.GLExplode(V,FVs,1.5,1.5,1.5,99,1));
GL.VIEW(GL.GLExplode(V,EVs,1.5,1.5,1.5,99,1));
meshes = GL.GLExplode(V,CVs[1:end],8,4,6,99,1);
GL.VIEW( push!( meshes, GL.GLFrame) );


text:
function ConstantRealDistribution(arg0::jdouble)
    return ConstantRealDistribution((jdouble,), arg0)
end

function cumulative_probability(obj::ConstantRealDistribution, arg0::jdouble)
    return jcall(obj, "cumulativeProbability", jdouble, (jdouble,), arg0)
end

function density(obj::ConstantRealDistribution, arg0::jdouble)
    return jcall(obj, "density", jdouble, (jdouble,), arg0)
end

function get_numerical_mean(obj::ConstantRealDistribution)
    return jcall(obj, "getNumericalMean", jdouble, ())
end

function get_numerical_variance(obj::ConstantRealDistribution)
    return jcall(obj, "getNumericalVariance", jdouble, ())
end

function get_support_lower_bound(obj::ConstantRealDistribution)
    return jcall(obj, "getSupportLowerBound", jdouble, ())
end

function get_support_upper_bound(obj::ConstantRealDistribution)
    return jcall(obj, "getSupportUpperBound", jdouble, ())
end

function inverse_cumulative_probability(obj::ConstantRealDistribution, arg0::jdouble)
    return jcall(obj, "inverseCumulativeProbability", jdouble, (jdouble,), arg0)
end

function is_support_connected(obj::ConstantRealDistribution)
    return jcall(obj, "isSupportConnected", jboolean, ())
end


text:
# ---
# title: 227. Basic Calculator II
# id: problem227
# author: Tian Jun
# date: 2020-10-31
# difficulty: Medium
# categories: String, Stack
# link: <https://leetcode.com/problems/basic-calculator-ii/description/>
# hidden: true
# ---
#
# Given a string `s` which represents an expression, _evaluate this expression
# and return its value_.
#
# The integer division should truncate toward zero.
#
#
#
# **Example 1:**
#
#
#
#     Input: s = "3+2*2"
#     Output: 7
#
#
# **Example 2:**
#
#
#
#     Input: s = " 3/2 "
#     Output: 1
#
#
# **Example 3:**
#
#
#
#     Input: s = " 3+5 / 2 "
#     Output: 5
#
#
#
#
# **Constraints:**
#
#   * `1 <= s.length <= 3 * 105`
#   * `s` consists of integers and operators `('+', '-', '*', '/')` separated by some number of spaces.
#   * `s` represents **a valid expression**.
#   * All the integers in the expression are non-negative integers in the range `[0, 231 - 1]`.
#   * The answer is **guaranteed** to fit in a **32-bit integer**.
#
#
## @lc code=start
using LeetCode

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


text:
"""Create the covariance matrix between the matrix X1 and X2 with the covariance function `kernel`"""
function kernelmatrix(X1::AbstractArray{T1},X2::AbstractArray{T2},kernel::Kernel) where {T1<:Real,T2<:Real}
    K = pairwise(getmetric(kernel),X1,X2,dims=1)
    v = getvariance(kernel)
    return lmul!(v,map!(kappa(kernel),K,K))
end

"""Compute the covariance matrix between the matrix X1 and X2 with the covariance function `kernel` in preallocated matrix K"""
function kernelmatrix!(K::AbstractArray{T1},X1::AbstractArray{T2},X2::AbstractArray{T3},kernel::Kernel) where {T1<:Real,T2<:Real,T3<:Real}
    (n1,n2) = size(K)
    @assert n1==size(X1,1)
    @assert n2==size(X2,1)
    pairwise!(K,getmetric(kernel),X1,X2,dims=1)
    v = getvariance(kernel)
    map!(kappa(kernel),K,K)
    return lmul!(v,K)
end

"""Compute the covariance matrix of the matrix X, optionally only compute the diagonal terms"""
function kernelmatrix(X::AbstractArray{T1},kernel::Kernel;diag::Bool=false) where {T1<:Real}
    if diag
        return kerneldiagmatrix(X,kernel)
    end
    K = pairwise(getmetric(kernel),X,dims=1)
    v = getvariance(kernel)
    # return v.*map(kappa(kernel),K)
    return v*map!(kappa(kernel),K,K)
end

"""Compute the covariance matrix of the matrix X in preallocated matrix K, optionally only compute the diagonal terms"""
function kernelmatrix!(K::AbstractArray{T1},X::AbstractArray{T2},kernel::Kernel; diag::Bool=false) where {T1<:Real,T2<:Real}
    if diag
        kerneldiagmatrix!(K,X,kernel)
    end
    (n1,n2) = size(K)
    @assert n1==size(X,1)
    @assert n1==n2
    pairwise!(K,getmetric(kernel),X,dims=1)
    v = getvariance(kernel)
    map!(kappa(kernel),K,K)
    return lmul!(v,K)
end

"""Compute only the diagonal elements of the covariance matrix"""
function kerneldiagmatrix(X::AbstractArray{T1},kernel::Kernel) where {T1<:Real}
    n = size(X,1)
    K = zeros(T1,n)
    v = getvariance(kernel)
    f = kappa(kernel)
    @simd for i in eachindex(K)
        @inbounds K[i] = f(evaluate(getmetric(kernel),X[i,:],X[i,:]))
    end
    return lmul!(v,K)
end

"""Compute only the diagonal elements of the covariance matrix in preallocated vector K"""
function kerneldiagmatrix!(K::AbstractVector{T1},X::AbstractArray{T2},kernel::Kernel) where {T1<:Real,T2<:Real}
    @assert size(K,1) == size(X,1)
    v = getvariance(kernel)
    f = kappa(kernel)
    @simd for i in eachindex(K)
        @inbounds K[i] = f(evaluate(getmetric(kernel),X[i,:],X[i,:]))
    end
    return lmul!(v,K)
end

"""Remapping of the gradients into a matrix with only 1 column and 1 row being non-zero"""
function CreateColumnRowMatrix(n::Int,iter::Int,gradient::AbstractVector{T1}) where {T1<:Real}
    K = zeros(T1,n,n)
    K[iter,:] .= gradient; K[:,iter] .= gradient;
    return K
end

"""Remapping of the gradients into a matrix with only 1 column being non-zero"""
function CreateColumnMatrix(n::Int,m::Int,iter::Int,gradient::AbstractVector{T1}) where {T1 <:Real}
    K = zeros(T1,n,m)
    K[:,iter] .= gradient;
    return K
end

"Compute the gradients given the inducing point locations, (general gradients are computed to be then remapped correctly)"
function computeIndPointsJ(model,iter::Int)
    Dnm = computeIndPointsJnm(model.kernel,model.X[model.MBIndices,:],model.inducingPoints[iter,:],iter,model.Knm)
    Dmm = computeIndPointsJmm(model.kernel,model.inducingPoints,iter,model.Kmm)
    Jnm = zeros(model.nDim,model.nSamplesUsed,model.m)
    Jmm = zeros(model.nDim,model.m,model.m)
    @inbounds for i in 1:model.nDim
        Jnm[i,:,:] .= CreateColumnMatrix(model.nSamplesUsed,model.m,iter,Dnm[:,i])
        Jmm[i,:,:] .= CreateColumnRowMatrix(model.m,iter,Dmm[:,i])
    end
    return Jnm,Jmm
    #Return dim*K*K tensors for computing the gradient
end