text:
# To run this script, activate the environment in tutorials:
# ]activate .

using Random, Statistics, LaTeXStrings, DataFrames, BenchmarkTools, CairoMakie, CSV
using Cropbox, LeafGasExchange, RCall, PlantBiophysics

# Parameters :
Random.seed!(1)
microbenchmark_steps = 100 # Number of times the microbenchmark is run
microbenchmark_evals = 1 # Number of times each sample is run to be sure of the output
N = 100        # Number of timesteps simulated for each microbenchmark step

# Create the ranges of input parameters
length_range = 10000
Rs = range(10, 500, length=length_range)
Ta = range(18, 40, length=length_range)
Wind = range(0.5, 20, length=length_range)
P = range(90, 101, length=length_range)
Rh = range(0.1, 0.98, length=length_range)
Ca = range(360, 900, length=length_range)
skyF = range(0.0, 1.0, length=length_range)
d = range(0.001, 0.5, length=length_range)
Jmax = range(200.0, 300.0, length=length_range)
Vmax = range(150.0, 250.0, length=length_range)
Rd = range(0.3, 2.0, length=length_range)
TPU = range(5.0, 20.0, length=length_range)
g0 = range(0.001, 2.0, length=length_range)
g1 = range(0.5, 15.0, length=length_range)
vars = hcat([Ta, Wind, P, Rh, Ca, Jmax, Vmax, Rd, Rs, skyF, d, TPU, g0, g1])

# Create the randomly drawn input set
set = [rand.(vars) for i in 1:N]
set = reshape(vcat(set...), (length(set[1]), length(set)))'
name = ["T", "Wind", "P", "Rh", "Ca", "JMaxRef", "VcMaxRef", "RdRef", "Rs", "sky_fraction", "d", "TPURef", "g0", "g1"]
set = DataFrame(set, name)
@. set[!, :vpd] = e_sat(set.T) - vapor_pressure(set.T, set.Rh)
@. set[!, :PPFD] = set.Rs * 0.48 * 4.57

######################################################################################################
# BENCHMARKING PLANTECOPHYS
######################################################################################################

# Import dependencies, install if not installed:
R"""
if(!require("plantecophys")){
    install.packages("plantecophys", repos = "https://cloud.r-project.org")
}
if(!require("microbenchmark")){
    install.packages("microbenchmark", repos = "https://cloud.r-project.org")
}
"""

# Make variables available to the R session
@rput set N microbenchmark_steps

# Benchmark:
R"""
# Define the function call in a function that takes a list as input to limit DataFrame overhead
function_EB <- function(input) {
    PhotosynEB(
        Tair = input$Tair, VPD = input$VPD, Wind = input$Wind,
        Wleaf = input$Wleaf,Ca = input$Ca,  StomatalRatio = 1,
        LeafAbs = input$LeafAbs, gsmodel = "BBOpti", g0 = input$g0, g1 = input$g1,
        alpha = 0.24, theta = 0.7, Jmax = input$Jmax,
        Vcmax = input$Vcmax, TPU = input$TPU, Rd = input$Rd,
        RH = input$RH, PPFD=input$PPFD, Patm = input$Patm
    )
}

time_PE = c()
for(i in seq_len(N)){
    # Put the inputs into a vector to limit dataframe overhead:
    input = list(
        Tair = set$T[i], VPD = set$vpd[i], Wind = set$Wind[i], Wleaf = set$d[i],
        Ca = set$Ca[i], LeafAbs = set$sky_fraction[i], g0 = set$g0[i], g1 = set$g1[i],
        Jmax = set$JMaxRef[i], Vcmax = set$VcMaxRef[i], TPU = set$TPURef[i],
        Rd = set$RdRef[i], RH = set$Rh[i]*100, PPFD=set$PPFD[i],Patm = set$P[i]
    )

    m = microbenchmark(function_EB(input), times = microbenchmark_steps)

    time_PE = append(time_PE,m$time * 10e-9) # transform in seconds
}
"""

@rget time_PE

######################################################################################################
# BENCHMARKING LEAFGASEXCHANGE.JL
######################################################################################################

time_LG = []
n_lg = fill(0, N)
for i in 1:N
    config = :Weather => (
        PFD=set.PPFD[i],
        CO2=set.Ca[i],
        RH=set.Rh[i] * 100,
        T_air=set.T[i],
        wind=set.Wind[i],
        P_air=set.P[i],
        g0=set.g0[i],
        g1=set.g1[i],
        Vcmax=set.VcMaxRef[i],
        Jmax=set.JMaxRef[i],
        Rd=set.RdRef[i],
        TPU=set.TPURef[i]
    )
    b_LG = @benchmark simulate($ModelC3MD; config=$config) evals = microbenchmark_evals samples = microbenchmark_steps
    append!(time_LG, b_LG.times .* 1e-9) # transform in seconds
    n_lg[i] = 1
end

######################################################################################################
# BENCHMARKING PLANTBIOPHYSICS.JL
######################################################################################################

time_PB = []
for i in 1:N
    leaf = LeafModels(
        energy=Monteith(),
        photosynthesis=Fvcb(
            VcMaxRef=set.VcMaxRef[i],
            JMaxRef=set.JMaxRef[i],
            RdRef=set.RdRef[i],
            TPURef=set.TPURef[i]
        ),
        stomatal_conductance=Medlyn(set.g0[i], set.g1[i]),
        Rₛ=set.Rs[i], sky_fraction=set.sky_fraction[i], PPFD=set.PPFD[i], d=set.d[i]
    )
    meteo = Atmosphere(T=set.T[i], Wind=set.Wind[i], P=set.P[i], Rh=set.Rh[i], Cₐ=set.Ca[i])

    b_PB = @benchmark energy_balance!($leaf, $meteo) evals = microbenchmark_evals samples = microbenchmark_steps
    append!(time_PB, b_PB.times .* 1e-9) # transform in seconds
end

######################################################################################################
# STATISTICS
######################################################################################################

# Defining a structure to hold basic statistics
struct StatResults
    mean::AbstractFloat
    median::AbstractFloat
    stddev::AbstractFloat
    min::AbstractFloat
    max::AbstractFloat
end

function Base.show(io::IO, ::MIME"text/plain", m::StatResults)
    print(
        io,
        "Benchmark:",
        "\nMean time -> ", m.mean, " ± ", m.stddev,
        "\nMedian time -> ", m.median,
        "\nMinimum time -> ", m.min,
        "\nMaximum time -> ", m.max,
    )
end
Base.show(io::IO, m::StatResults) = print(io, m.mean, "(±", m.stddev, ')')

# Function computing the basic statistics
function basic_stat(df)
    m = mean(df)
    med = median(df)
    std = Statistics.std(df)
    min = findmin(df)[1]
    max = findmax(df)[1]
    return StatResults(m, med, std, min, max)
end

statsPB = basic_stat(time_PB)
statsPE = basic_stat(time_PE)
statsLG = basic_stat(time_LG)

factorPE = mean(time_PE) / mean(time_PB)
factorLG = mean(time_LG) / mean(time_PB)

######################################################################################################
# PLOTTING
######################################################################################################

function plot_benchmark_Makie(statsPB, statsPE, statsLG, time_PB, time_PE, time_LG)
    size_inches = (6.7, 5)
    size_pt = 72 .* size_inches
    bins = 220
    noto_sans = assetpath("fonts", "NotoSans-Regular.ttf")
    fig = Figure(backgroundcolor=RGBf(1, 1, 1), resolution=size_pt, font=noto_sans, fontsize=10)
    ep = 1e-9
    extr = extrema(vcat(time_PB, time_PE, time_LG))
    interval = (extr[1] * 1e-1, extr[2])

    ga = fig[1, 1] = GridLayout()

    axa = Axis(ga[1, 1], title="(a) PlantBiophysics.jl", xscale=log10, titlealign=:left, titlesize=10)
    stddevi = poly!(
        axa,
        Rect(max(ep, statsPB.mean - statsPB.stddev), 0.0, 2 * statsPB.stddev, 1),
        color=(:orange, 0.3), yautolimits=false
    )
    moy = vlines!(axa, statsPB.mean; color=:red, linewidth=3, linestyle=:dot)
    hist!(axa, time_PB, normalization=:probability, bins=bins)
    # h = axa.finallimits[].widths[2]
    axislegend(axa, [stddevi, moy], ["95% confidence interval", "Mean"], "", position=:rb,
        orientation=:vertical, labelsize=8, framevisible=false)
    xlims!(axa, interval)

    axb = Axis(ga[2, 1], title="(b) plantecophys", xscale=log10, ylabel="Density", titlealign=:left, titlesize=10)
    stddevi = poly!(
        axb,
        Rect(statsPE.mean - statsPE.stddev, 0.0, 2 * statsPE.stddev, 1),
        color=(:orange, 0.3), yautolimits=false
    )
    vlines!(axb, statsPE.mean; color=:red, linewidth=3, linestyle=:dot)
    hist!(axb, time_PE, normalization=:probability, bins=bins)
    xlims!(axb, interval)

    # axc = Axis(ga[3, 1], title="(c) LeafGasExchange.jl", yminorticks=IntervalsBetween(10),
    #     xscale=log10, xminorticks=IntervalsBetween(10), yminorgridvisible=true, yminorticksvisible=true,
    #     xminorgridvisible=true, xminorticksvisible=true, xlabel="Time (s)")
    axc = Axis(ga[3, 1], title="(c) LeafGasExchange.jl", xscale=log10, xlabel="Time (s)", titlealign=:left, titlesize=10)
    stddevi = poly!(
        axc,
        Rect(statsLG.mean - statsLG.stddev, 0.0, 2 * statsLG.stddev, 1),
        color=(:orange, 0.3), yautolimits=false
    )
    vlines!(axc, statsLG.mean; color=:red, linewidth=3, linestyle=:dot)
    hist!(axc, time_LG, normalization=:probability, bins=bins)
    xlims!(axc, interval)

    rowgap!(ga, 7)
    hidexdecorations!(axa, grid=false)
    hidexdecorations!(axb, grid=false)
    fig
end

######################################################################################################
# SAVING
######################################################################################################

# Save the figure:
fig = plot_benchmark_Makie(statsPB, statsPE, statsLG, time_PB, time_PE, time_LG)
save("benchmark_each_time_steps.png", fig, px_per_unit=3)

# Write overall timings:
df = DataFrame(
    [getfield(j, i) for i in fieldnames(StatResults), j in [statsPB, statsPE, statsLG]],
    ["PlantBiophysics", "plantecophys", "LeafGasExchange"]
)
insertcols!(df, 1, :Stat => [fieldnames(StatResults)...])
CSV.write("benchmark.csv", df)

# Write timing for each sample:
CSV.write("benchmark_full.csv",
    DataFrame(
        "package" => vcat(
            [
                repeat([i.first], length(i.second)) for i in [
                    "PlantBiophysics" => time_PB,
                    "plantecophys" => time_PE,
                    "LeafGasExchange" => time_LG
                ]
            ]...
        ),
        "sample_time" => vcat(time_PB, time_PE, time_LG)
    )
)

######################################################################################################
# Benchmark with limiting overhead
######################################################################################################

# Benchmark but passing all time-steps at once and computing the average simulation time.
# This is helpfull to avoid the overhead of calling the wrapper function, the downside is that
# we don't have the variability between simulations

# Parameters :
# Random.seed!(1)
# microbenchmark_steps2 = 1000 # Number of times the microbenchmark is run
# microbenchmark_evals2 = 1 # Number of times each sample is run to be sure of the output
# N2 = 100        # Number of timesteps simulated for each microbenchmark step

# # Create the randomly drawn input set
# set2 = [rand.(vars) for i in 1:N2]
# set2 = reshape(vcat(set2...), (length(set2[1]), length(set2)))'
# name = ["T", "Wind", "P", "Rh", "Ca", "JMaxRef", "VcMaxRef", "RdRef", "Rs", "sky_fraction", "d", "TPURef", "g0", "g1"]
# set2 = DataFrame(set2, name)
# @. set2[!, :vpd] = e_sat(set2.T) - vapor_pressure(set2.T, set2.Rh)
# @. set2[!, :PPFD] = set2.Rs * 0.48 * 4.57

# # Make variables available to the R session
# @rput set2 N2 microbenchmark_steps2
# R"""
# m2 = microbenchmark(
#     PhotosynEB(
#         Tair = set2$T, VPD = set2$vpd, Wind = set2$Wind,
#         Wleaf = set2$d, Ca = set2$Ca,  StomatalRatio = 1,
#         LeafAbs = set2$sky_fraction, gsmodel = "BBOpti", g0 = set2$g0, g1 = set2$g1,
#         alpha = 0.24, theta = 0.7, Jmax = set2$JMaxRef,
#         Vcmax = set2$VcMaxRef, TPU = set2$TPURef, Rd = set2$RdRef,
#         RH = set2$Rh, PPFD=set2$PPFD, Patm = set2$P
#     ),
#     times = microbenchmark_steps2
# )

#     time_PE2 = m2$time / N2 * 10e-9 # transform in seconds
# """
# @rget time_PE2

# # Configuration:
# leaves = [
#     begin
#         LeafModels(
#             energy=Monteith(),
#             photosynthesis=Fvcb(VcMaxRef=set2.VcMaxRef[i], JMaxRef=set2.JMaxRef[i], RdRef=set2.RdRef[i], TPURef=set2.TPURef[i]),
#             stomatal_conductance=Medlyn(set2.g0[i], set2.g1[i]),
#             Rₛ=set2.Rs[i], sky_fraction=set2.sky_fraction[i],
#             PPFD=set2.PPFD[i], d=set2.d[i]
#         )
#     end for i in 1:size(set2, 1)
# ]

# meteos = [
#     begin
#         Atmosphere(T=set2.T[i], Wind=set2.Wind[i], P=set2.P[i], Rh=set2.Rh[i], Cₐ=set2.Ca[i])
#     end for i in 1:size(set2, 1)
# ]

# function eval_PB(leaves, meteos)
#     for i in 1:N2
#         energy_balance!(leaves[i], meteos[i])
#     end
# end

# b_PB2 = @benchmark eval_PB($leaves, $meteos) evals = microbenchmark_evals2 samples = microbenchmark_steps2
# time_PB2 = b_PB2.times ./ N2 .* 1e-9 # transform in seconds

# # For LeafGasExchange:
# # This is usefull to know the overhead cost from Julia to Python and back.
# configs = [
#     begin
#         :Weather => (
#             PFD=set2.PPFD[i],
#             CO2=set2.Ca[i],
#             RH=set2.Rh[i] * 100,
#             T_air=set2.T[i],
#             wind=set2.Wind[i],
#             P_air=set2.P[i],
#             g0=set2.g0[i],
#             g1=set2.g1[i],
#             Vcmax=set2.VcMaxRef[i],
#             Jmax=set2.JMaxRef[i],
#             Rd=set2.RdRef[i],
#             TPU=set2.TPURef[i]
#         )
#     end for i in 1:size(set2, 1)
# ]
# b_LG2 = @benchmark simulate($ModelC3MD; configs=$configs) evals = microbenchmark_evals2 samples = microbenchmark_steps2
# time_LG2 = b_LG2.times ./ N2 .* 1e-9 # transform in seconds

# ######################################################################################################
# # SAVING
# ######################################################################################################

# statsPB2 = basic_stat(time_PB2)
# statsPE2 = basic_stat(time_PE2)
# statsLG2 = basic_stat(time_LG2)

# plot_benchmark(statsPB2, statsPE2, statsLG2)
# savefig("benchmark_all_steps.png")

max_stars_repo_name: VEZY/PlantBiophysics-paper

max_issues_repo_name: VEZY/PlantBiophysics-paper

max_forks_repo_name: VEZY/PlantBiophysics-paper


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

module Bar
import Baz, Foo
this = "Bar"
end

ext: jl

max_stars_repo_name: greimel/julia

max_issues_repo_name: greimel/julia

max_forks_repo_name: greimel/julia


text:
"""
Bridges `CP.BinPacking{CP.FIXED_CAPACITY_BINPACKING, T}` to
`CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING, T}` by creating capacity
variables.
"""
struct FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T} <: MOIBC.AbstractBridge
    capa_var::Vector{MOI.VariableIndex}
    capa_con::Vector{MOI.ConstraintIndex{MOI.VariableIndex, MOI.Integer}}
    capa_bound::Vector{MOI.ConstraintIndex{MOI.ScalarAffineFunction{T}, MOI.LessThan{T}}}
    bp::MOI.ConstraintIndex{MOI.VectorAffineFunction{T}, CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING, T}}
end

function MOIBC.bridge_constraint(
    ::Type{FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T}},
    model,
    f::MOI.VectorOfVariables,
    s::CP.BinPacking{CP.FIXED_CAPACITY_BINPACKING, T},
) where {T}
    return MOIBC.bridge_constraint(
        FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
        model,
        MOI.VectorAffineFunction{T}(f),
        s,
    )
end

function MOIBC.bridge_constraint(
    ::Type{FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T}},
    model,
    f::MOI.VectorAffineFunction{T},
    s::CP.BinPacking{CP.FIXED_CAPACITY_BINPACKING, T},
) where {T <: Integer}
    # Add the capacity integer variables (because its values are ensured to
    # be integers, T <: Integer).
    capa_var, capa_con = MOI.add_constrained_variables(model, [MOI.Integer() for _ in 1:s.n_bins])

    # Add the capacity constraints.
    capa_bound = [
        MOI.add_constraint(
            model,
            one(T) * capa_var[bin],
            MOI.LessThan(s.capacities[bin])
        )
        for bin in 1:s.n_bins
    ]

    # Create the variable-capacity bin-packing constraint.
    f_scalars = MOIU.scalarize(f)
    new_f = MOIU.vectorize(
        MOI.ScalarAffineFunction{T}[
            f_scalars[1:s.n_bins]...,
            capa_var...,
            f_scalars[s.n_bins+1:end]...
        ]
    )
    bp_set = CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING}(s.n_bins, s.n_items, s.weights)
    bp = MOI.add_constraint(model, new_f, bp_set)

    return FixedCapacityBinPacking2VariableCapacityBinPackingBridge(capa_var, capa_con, capa_bound, bp)
end

function MOIBC.bridge_constraint(
    ::Type{FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T}},
    model,
    f::MOI.VectorAffineFunction{T},
    s::CP.BinPacking{CP.FIXED_CAPACITY_BINPACKING, T},
) where {T <: Real}
    # Add the capacity variables, without integer constraint.
    capa_var = MOI.add_variables(model, s.n_bins)
    capa_con = MOI.ConstraintIndex{MOI.VariableIndex, MOI.Integer}[]

    # Add the capacity constraints.
    capa_bound = [
        MOI.add_constraint(
            model,
            one(T) * capa_var[bin],
            MOI.LessThan(s.capacities[bin])
        )
        for bin in 1:s.n_bins
    ]

    # Create the variable-capacity bin-packing constraint.
    f_scalars = MOIU.scalarize(f)
    new_f = MOIU.vectorize(
        MOI.ScalarAffineFunction{T}[
            f_scalars[1:s.n_bins]...,
            capa_var...,
            f_scalars[s.n_bins+1:end]...
        ]
    )
    bp_set = CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING}(s.n_bins, s.n_items, s.weights)
    bp = MOI.add_constraint(model, new_f, bp_set)

    return FixedCapacityBinPacking2VariableCapacityBinPackingBridge(capa_var, capa_con, capa_bound, bp)
end

function MOI.supports_constraint(
    ::Type{FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T}},
    ::Union{Type{MOI.VectorOfVariables}, Type{MOI.VectorAffineFunction{T}}},
    ::Type{CP.BinPacking{CP.FIXED_CAPACITY_BINPACKING, T}},
) where {T}
    return true
end

function MOIB.added_constrained_variable_types(::Type{<:FixedCapacityBinPacking2VariableCapacityBinPackingBridge{<:Integer}})
    return [(MOI.Integer,)]
end

function MOIB.added_constrained_variable_types(::Type{<:FixedCapacityBinPacking2VariableCapacityBinPackingBridge{<:Real}})
    return Tuple{Type}[]
end

function MOIB.added_constraint_types(::Type{FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T}}) where {T}
    return [
        (MOI.VectorAffineFunction{T}, CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING, T}),
        (MOI.ScalarAffineFunction{T}, MOI.LessThan{T}),
    ]
end

function MOI.get(b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge, ::MOI.NumberOfVariables)
    return length(b.capa_var)
end

function MOI.get(
    ::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.NumberOfConstraints{
        MOI.VectorAffineFunction{T},
        CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING, T},
    },
) where {T}
    return 1
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.NumberOfConstraints{
        MOI.VariableIndex,
        MOI.Integer,
    },
) where {T}
    return length(b.capa_con)
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.NumberOfConstraints{
        MOI.ScalarAffineFunction{T},
        MOI.LessThan{T},
    },
) where {T}
    return length(b.capa_bound)
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.ListOfVariableIndices,
) where {T}
    return copy(b.capa_var)
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.ListOfConstraintIndices{
        MOI.VectorAffineFunction{T},
        CP.BinPacking{CP.VARIABLE_CAPACITY_BINPACKING, T},
    },
) where {T}
    return [b.bp]
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.ListOfConstraintIndices{
        MOI.VariableIndex,
        MOI.Integer,
    },
) where {T}
    return copy(b.capa_con)
end

function MOI.get(
    b::FixedCapacityBinPacking2VariableCapacityBinPackingBridge{T},
    ::MOI.ListOfConstraintIndices{
        MOI.ScalarAffineFunction{T},
        MOI.LessThan{T},
    },
) where {T}
    return copy(b.capa_bound)
end


text:
 #=
Author:= Shreeda Bhat
License:= MIT
UTF-8 =#

#=
In Julia the for loop is simple, Unlike python, matlab it won't use vectorized code.

    `for.....end`

=#

println("***********************1st Example*************************")
for num = 1:10 # here the colon is used represent the range like in python
    println("number is $num")
end

println("***********************2nd Example*************************")

x = ["shreeda",1,2,0.3,4]

for i in x
    println("what i'm passing now is $i")
end

#=
Note: in the condition it should be `variable = iterable`
=#

println("***********************3rd Example*************************")

for i in 1:10
    println(i)
end

println("***********************4th Example*************************")

my_bio = ["shreeda","20","Astaqc Consulting"]

for bio in my_bio
    println("let's see this array $bio")
end

println("*********************** 5th Example*************************")
m,n = 3,3
array = zeros(m,n)
println(array)

for i in 1:m
    for j in 1:n
        array[i,j]=i+j
    end
end
println(array)

println("*********************** 6th Example*************************")
# Let's do list comprehension
B = zeros(3,3)
for i in 1:m,j in 1:n
    B[i,j]=i+j
end
println(B)

println("*********************** 7th Example*************************")
# Sorry its array comprehension
c= [i+j for i in 1:n, j in 1:m]
print(c)


text:
using ImageRegistration
using Test

include("test_boundingbox.jl")
include("test_imwarp.jl")
include("test_incidence_to_triangles.jl")
include("test_mesh.jl")
include("test_meshwarp.jl")
include("test_registration.jl")
include("test_transforms.jl")
include("test_geometry.jl")


text:
extract_ps(root) = extract_ps!([], root)
extract_ps(doc::HTMLDocument) = extract_ps(doc.root[2])

extract_ps!(accum, elem::HTMLElement{:p}) = push!(accum, elem)
extract_ps!(accum, elem::HTMLElement) = (foreach(x -> extract_ps!(accum, x), elem.children); accum)

function parse(doc; detect_labels=true)
    stack = []
    questions = Question[]

    ps = extract_ps(doc)
    while !is_question(first(ps))
        popfirst!(ps)
    end

    question_counter=0
    response_counter=0

    question_labels = Set([])

    for item in ps
        if is_response(item)
            if !isempty(stack)
                push!(
                    questions,
                    question(
                        copy(stack),
                        [],
                        "$(length(questions)+1)",
                        detect_label = detect_labels
                    )
                )
                empty!(stack)
                response_counter = 'A'

                lab = last(questions).label
                lab ∈ question_labels && @warn "Duplicate label detected: $lab"
                push!(question_labels, lab)
            end

            r = response(item, string(response_counter))
            response_counter += 1

            push!(last(questions).responses, r)
        else
            push!(stack, item)
        end
    end

    return questions

end

read_key(fn) = OrderedDict(map(x -> Pair(split(x)...), readlines(fn)))

rekey(r::Response, correct_lab) =
    Response(r.content, r.label, r.id, r.label == correct_lab)
function rekey!(q::Question, key::AbstractDict)
    if haskey(key, q.label)
        correct = key[q.label]
        any(isequal(correct), getfield.(q.responses, :label)) ||
            @warn "Response not found: $correct"
        q.responses .= rekey.(q.responses, correct)
    else
        @warn "No entry found for question $(q.label)"
    end
end
rekey!(qs::Vector{Question}, key::AbstractDict) = (rekey!.(qs, Ref(key)); qs)

function cleanup_noneall(q::Question)
    has_none = any(occursin.(r"none of the above"i, q.responses))
    has_all = any(occursin.(r"all of the above"i, q.responses))

    if has_none && has_all
        for r in q.responses
            replace!(r, r"(none of) the above"i => s"\1 these")
            replace!(r, r"(all of) the above"i => s"\1 these (except \"none\")")
        end
    elseif has_none || has_all
        for r in q.responses
            replace!(r, r"(none of) the above"i => s"\1 these")
            replace!(r, r"(all of) the above"i => s"\1 these")
        end
    end

    return q
end


text:
module DesignAutomation

include("bracket.jl")
include("golden_section.jl")
include("ccs.jl")
include("hooke_jeeves.jl")
include("uninformed_search.jl")
include("airport.jl")
include("exhaustive_search.jl")
include("random_hill.jl")

export bracket, golden_section,
    ccs, hooke_jeeves,
    dfs, bfs,
    exhaustive_search, random_hill

end # module

max_stars_repo_name: camirmas/DesignAutomation

max_issues_repo_name: camirmas/DesignAutomation

max_forks_repo_name: camirmas/DesignAutomation


text:
using Test
using AbstractTrees

using FlatZincParser
using FlatZincParser: TokenStream, AbstractNode, DataNode, Node, next!, reset!, context
using FlatZincParser: match!, match_many!, match_token, match_many

include("lexer_test.jl")
include("parser_test.jl")


text:
# this file not included in package

struct PSF_gauss2D <: PSF
    σ::Float32
end

function genmodel!(s::StateFlatBg, sz::Int32, psf::PSF_gauss2D, model::Array{Float32,2})
    for ii = 1:sz
        for jj = 1:sz
            model[ii,jj] = s.bg
            for nn = 1:s.n
                model[ii,jj] += s.photons[nn] / (2 * π * psf.σ^2) *
                exp(-(ii - s.y[nn])^2 / (2 * psf.σ^2)) *
                exp(-(jj - s.x[nn])^2 / (2 * psf.σ^2))
            end
        end
    end
end

function genmodel_gauss2D_CUDA!(s_n::Int32,s_x, s_y,
    s_photons,s_bg::Float32, sz::Int32, σ::Float32, model)
     #note that the 2d array is linearized and using 1-based indexing in kernel
     ii = blockIdx().x
     jj = threadIdx().x

     idx=(ii-1)*sz+jj
     model[idx] = s_bg + 1f-4
     for nn = 1:s_n
         model[idx] += s_photons[nn] / (2 * π * σ^2) *
                 exp(-(ii - s_y[nn])^2 / (2 * σ^2)) *
                 exp(-(jj - s_x[nn])^2 / (2 * σ^2))
     end
     return nothing
end

function genmodel!(s::StateFlatBg, sz::Int32, psf::PSF_gauss2D, model::CuArray{Float32,2})
    @cuda threads=sz blocks=sz genmodel_gauss2D_CUDA!(s.n,s.x,s.y,s.photons,s.bg, sz, psf.σ, model)
end

function max2int(psf::PSF_gauss2D)
    return 2*pi*psf.σ
end


text:
# This file is a part of Julia. License is MIT: https://julialang.org/license

# Document NTuple here where we have everything needed for the doc system
"""
    NTuple{N, T}

A compact way of representing the type for a tuple of length `N` where all elements are of type `T`.

# Examples
```jldoctest
julia> isa((1, 2, 3, 4, 5, 6), NTuple{6, Int})
true
```
"""
NTuple

## indexing ##

length(@nospecialize t::Tuple) = nfields(t)
firstindex(@nospecialize t::Tuple) = 1
lastindex(@nospecialize t::Tuple) = length(t)
size(@nospecialize(t::Tuple), d::Integer) = (d == 1) ? length(t) : throw(ArgumentError("invalid tuple dimension $d"))
axes(@nospecialize t::Tuple) = (OneTo(length(t)),)
@eval getindex(@nospecialize(t::Tuple), i::Int) = getfield(t, i, $(Expr(:boundscheck)))
@eval getindex(@nospecialize(t::Tuple), i::Real) = getfield(t, convert(Int, i), $(Expr(:boundscheck)))
getindex(t::Tuple, r::AbstractArray{<:Any,1}) = ([t[ri] for ri in r]...,)
getindex(t::Tuple, b::AbstractArray{Bool,1}) = length(b) == length(t) ? getindex(t, findall(b)) : throw(BoundsError(t, b))
getindex(t::Tuple, c::Colon) = t

# returns new tuple; N.B.: becomes no-op if i is out-of-bounds
setindex(x::Tuple, v, i::Integer) = (@_inline_meta; _setindex(v, i, x...))
function _setindex(v, i::Integer, first, tail...)
    @_inline_meta
    return (ifelse(i == 1, v, first), _setindex(v, i - 1, tail...)...)
end
_setindex(v, i::Integer) = ()

## iterating ##

function iterate(@nospecialize(t::Tuple), i::Int=1)
    @_inline_meta
    return (1 <= i <= length(t)) ? (@inbounds t[i], i + 1) : nothing
end

keys(@nospecialize t::Tuple) = OneTo(length(t))

prevind(@nospecialize(t::Tuple), i::Integer) = Int(i)-1
nextind(@nospecialize(t::Tuple), i::Integer) = Int(i)+1

function keys(t::Tuple, t2::Tuple...)
    @_inline_meta
    OneTo(_maxlength(t, t2...))
end
_maxlength(t::Tuple) = length(t)
function _maxlength(t::Tuple, t2::Tuple, t3::Tuple...)
    @_inline_meta
    max(length(t), _maxlength(t2, t3...))
end

# this allows partial evaluation of bounded sequences of next() calls on tuples,
# while reducing to plain next() for arbitrary iterables.
indexed_iterate(t::Tuple, i::Int, state=1) = (@_inline_meta; (getfield(t, i), i+1))
indexed_iterate(a::Array, i::Int, state=1) = (@_inline_meta; (a[i], i+1))
function indexed_iterate(I, i)
    x = iterate(I)
    x === nothing && throw(BoundsError(I, i))
    x
end
function indexed_iterate(I, i, state)
    x = iterate(I, state)
    x === nothing && throw(BoundsError(I, i))
    x
end

# Use dispatch to avoid a branch in first
first(::Tuple{}) = throw(ArgumentError("tuple must be non-empty"))
first(t::Tuple) = t[1]

# eltype

eltype(::Type{Tuple{}}) = Bottom
function eltype(t::Type{<:Tuple{Vararg{E}}}) where {E}
    if @isdefined(E)
        return E
    else
        # TODO: need to guard against E being miscomputed by subtyping (ref #23017)
        # and compute the result manually in this case
        return _compute_eltype(t)
    end
end
eltype(t::Type{<:Tuple}) = _compute_eltype(t)
function _compute_eltype(t::Type{<:Tuple})
    @_pure_meta
    t isa Union && return promote_typejoin(eltype(t.a), eltype(t.b))
    t´ = unwrap_unionall(t)
    r = Union{}
    for ti in t´.parameters
        r = promote_typejoin(r, rewrap_unionall(unwrapva(ti), t))
    end
    return r
end

# version of tail that doesn't throw on empty tuples (used in array indexing)
safe_tail(t::Tuple) = tail(t)
safe_tail(t::Tuple{}) = ()

# front (the converse of tail: it skips the last entry)

"""
    front(x::Tuple)::Tuple

Return a `Tuple` consisting of all but the last component of `x`.

# Examples
```jldoctest
julia> Base.front((1,2,3))
(1, 2)

julia> Base.front(())
ERROR: ArgumentError: Cannot call front on an empty tuple.
```
"""
function front(t::Tuple)
    @_inline_meta
    _front(t...)
end
_front() = throw(ArgumentError("Cannot call front on an empty tuple."))
_front(v) = ()
function _front(v, t...)
    @_inline_meta
    (v, _front(t...)...)
end

## mapping ##

# 1 argument function
map(f, t::Tuple{})              = ()
map(f, t::Tuple{Any,})          = (f(t[1]),)
map(f, t::Tuple{Any, Any})      = (f(t[1]), f(t[2]))
map(f, t::Tuple{Any, Any, Any}) = (f(t[1]), f(t[2]), f(t[3]))
map(f, t::Tuple)                = (@_inline_meta; (f(t[1]), map(f,tail(t))...))
# stop inlining after some number of arguments to avoid code blowup
const Any16{N} = Tuple{Any,Any,Any,Any,Any,Any,Any,Any,
                       Any,Any,Any,Any,Any,Any,Any,Any,Vararg{Any,N}}
const All16{T,N} = Tuple{T,T,T,T,T,T,T,T,
                         T,T,T,T,T,T,T,T,Vararg{T,N}}
function map(f, t::Any16)
    n = length(t)
    A = Vector{Any}(undef, n)
    for i=1:n
        A[i] = f(t[i])
    end
    (A...,)
end
# 2 argument function
map(f, t::Tuple{},        s::Tuple{})        = ()
map(f, t::Tuple{Any,},    s::Tuple{Any,})    = (f(t[1],s[1]),)
map(f, t::Tuple{Any,Any}, s::Tuple{Any,Any}) = (f(t[1],s[1]), f(t[2],s[2]))
function map(f, t::Tuple, s::Tuple)
    @_inline_meta
    (f(t[1],s[1]), map(f, tail(t), tail(s))...)
end
function map(f, t::Any16, s::Any16)
    n = length(t)
    A = Vector{Any}(undef, n)
    for i = 1:n
        A[i] = f(t[i], s[i])
    end
    (A...,)
end
# n argument function
heads(ts::Tuple...) = map(t -> t[1], ts)
tails(ts::Tuple...) = map(tail, ts)
map(f, ::Tuple{}...) = ()
function map(f, t1::Tuple, t2::Tuple, ts::Tuple...)
    @_inline_meta
    (f(heads(t1, t2, ts...)...), map(f, tails(t1, t2, ts...)...)...)
end
function map(f, t1::Any16, t2::Any16, ts::Any16...)
    n = length(t1)
    A = Vector{Any}(undef, n)
    for i = 1:n
        A[i] = f(t1[i], t2[i], map(t -> t[i], ts)...)
    end
    (A...,)
end

# mapafoldl, based on afold in operators.jl
mapafoldl(F,op,a) = a
mapafoldl(F,op,a,b) = op(a,F(b))
mapafoldl(F,op,a,b,c...) = mapafoldl(F, op, op(a,F(b)), c...)
function mapafoldl(F,op,a,b,c,d,e,f,g,h,i,j,k,l,m,n,o,p,qs...)
    y = op(op(op(op(op(op(op(op(op(op(op(op(op(op(op(a,F(b)),F(c)),F(d)),F(e)),F(f)),F(g)),F(h)),F(i)),F(j)),F(k)),F(l)),F(m)),F(n)),F(o)),F(p))
    for x in qs; y = op(y,F(x)); end
    y
end
mapfoldl_impl(f, op, nt::NamedTuple{(:init,)}, t::Tuple) = mapafoldl(f, op, nt.init, t...)
mapfoldl_impl(f, op, nt::NamedTuple{()}, t::Tuple) = mapafoldl(f, op, f(t[1]), tail(t)...)
mapfoldl_impl(f, op, nt::NamedTuple{()}, t::Tuple{}) = mapreduce_empty_iter(f, op, t, IteratorEltype(t))

# type-stable padding
fill_to_length(t::NTuple{N,Any}, val, ::Val{N}) where {N} = t
fill_to_length(t::Tuple{}, val, ::Val{1}) = (val,)
fill_to_length(t::Tuple{Any}, val, ::Val{2}) = (t..., val)
fill_to_length(t::Tuple{}, val, ::Val{2}) = (val, val)
#function fill_to_length(t::Tuple, val, ::Val{N}) where {N}
#    @_inline_meta
#    return (t..., ntuple(i -> val, N - length(t))...)
#end

# constructing from an iterator

# only define these in Base, to avoid overwriting the constructors
# NOTE: this means this constructor must be avoided in Core.Compiler!
if nameof(@__MODULE__) === :Base

(::Type{T})(x::Tuple) where {T<:Tuple} = convert(T, x)  # still use `convert` for tuples

(::Type{T})(itr) where {T<:Tuple} = _totuple(T, itr)

_totuple(::Type{Tuple{}}, itr, s...) = ()

function _totuple_err(@nospecialize T)
    @_noinline_meta
    throw(ArgumentError("too few elements for tuple type $T"))
end

function _totuple(T, itr, s...)
    @_inline_meta
    y = iterate(itr, s...)
    y === nothing && _totuple_err(T)
    (convert(tuple_type_head(T), y[1]), _totuple(tuple_type_tail(T), itr, y[2])...)
end

# use iterative algorithm for long tuples
function _totuple(T::Type{All16{E,N}}, itr) where {E,N}
    len = N+16
    elts = collect(E, Iterators.take(itr,len))
    if length(elts) != len
        _totuple_err(T)
    end
    (elts...,)
end

_totuple(::Type{Tuple{Vararg{E}}}, itr, s...) where {E} = (collect(E, Iterators.rest(itr,s...))...,)

_totuple(::Type{Tuple}, itr, s...) = (collect(Iterators.rest(itr,s...))...,)

end

## comparison ##

isequal(t1::Tuple, t2::Tuple) = (length(t1) == length(t2)) && _isequal(t1, t2)
_isequal(t1::Tuple{}, t2::Tuple{}) = true
_isequal(t1::Tuple{Any}, t2::Tuple{Any}) = isequal(t1[1], t2[1])
_isequal(t1::Tuple, t2::Tuple) = isequal(t1[1], t2[1]) && _isequal(tail(t1), tail(t2))
function _isequal(t1::Any16, t2::Any16)
    for i = 1:length(t1)
        if !isequal(t1[i], t2[i])
            return false
        end
    end
    return true
end

==(t1::Tuple, t2::Tuple) = (length(t1) == length(t2)) && _eq(t1, t2)
_eq(t1::Tuple{}, t2::Tuple{}) = true
_eq_missing(t1::Tuple{}, t2::Tuple{}) = missing
function _eq(t1::Tuple, t2::Tuple)
    eq = t1[1] == t2[1]
    if eq === false
        return false
    elseif ismissing(eq)
        return _eq_missing(tail(t1), tail(t2))
    else
        return _eq(tail(t1), tail(t2))
    end
end
function _eq_missing(t1::Tuple, t2::Tuple)
    eq = t1[1] == t2[1]
    if eq === false
        return false
    else
        return _eq_missing(tail(t1), tail(t2))
    end
end
function _eq(t1::Any16, t2::Any16)
    anymissing = false
    for i = 1:length(t1)
        eq = (t1[i] == t2[i])
        if ismissing(eq)
            anymissing = true
        elseif !eq
           return false
       end
    end
    return anymissing ? missing : true
end

const tuplehash_seed = UInt === UInt64 ? 0x77cfa1eef01bca90 : 0xf01bca90
hash(::Tuple{}, h::UInt) = h + tuplehash_seed
hash(t::Tuple, h::UInt) = hash(t[1], hash(tail(t), h))
function hash(t::Any16, h::UInt)
    out = h + tuplehash_seed
    for i = length(t):-1:1
        out = hash(t[i], out)
    end
    return out
end

<(::Tuple{}, ::Tuple{}) = false
<(::Tuple{}, ::Tuple) = true
<(::Tuple, ::Tuple{}) = false
function <(t1::Tuple, t2::Tuple)
    a, b = t1[1], t2[1]
    eq = (a == b)
    if ismissing(eq)
        return missing
    elseif !eq
        return a < b
    end
    return tail(t1) < tail(t2)
end
function <(t1::Any16, t2::Any16)
    n1, n2 = length(t1), length(t2)
    for i = 1:min(n1, n2)
        a, b = t1[i], t2[i]
        eq = (a == b)
        if ismissing(eq)
            return missing
        elseif !eq
           return a < b
        end
    end
    return n1 < n2
end

isless(::Tuple{}, ::Tuple{}) = false
isless(::Tuple{}, ::Tuple) = true
isless(::Tuple, ::Tuple{}) = false

"""
    isless(t1::Tuple, t2::Tuple)

Returns true when t1 is less than t2 in lexicographic order.
"""
function isless(t1::Tuple, t2::Tuple)
    a, b = t1[1], t2[1]
    isless(a, b) || (isequal(a, b) && isless(tail(t1), tail(t2)))
end
function isless(t1::Any16, t2::Any16)
    n1, n2 = length(t1), length(t2)
    for i = 1:min(n1, n2)
        a, b = t1[i], t2[i]
        if !isequal(a, b)
            return isless(a, b)
        end
    end
    return n1 < n2
end

## functions ##

isempty(x::Tuple{}) = true
isempty(@nospecialize x::Tuple) = false

revargs() = ()
revargs(x, r...) = (revargs(r...)..., x)

reverse(t::Tuple) = revargs(t...)

## specialized reduction ##

# TODO: these definitions cannot yet be combined, since +(x...)
# where x might be any tuple matches too many methods.
# TODO: this is inconsistent with the regular sum in cases where the arguments
# require size promotion to system size.
sum(x::Tuple{Any, Vararg{Any}}) = +(x...)

# NOTE: should remove, but often used on array sizes
# TODO: this is inconsistent with the regular prod in cases where the arguments
# require size promotion to system size.
prod(x::Tuple{}) = 1
prod(x::Tuple{Any, Vararg{Any}}) = *(x...)

all(x::Tuple{}) = true
all(x::Tuple{Bool}) = x[1]
all(x::Tuple{Bool, Bool}) = x[1]&x[2]
all(x::Tuple{Bool, Bool, Bool}) = x[1]&x[2]&x[3]
# use generic reductions for the rest

any(x::Tuple{}) = false
any(x::Tuple{Bool}) = x[1]
any(x::Tuple{Bool, Bool}) = x[1]|x[2]
any(x::Tuple{Bool, Bool, Bool}) = x[1]|x[2]|x[3]

# equivalent to any(f, t), to be used only in bootstrap
_tuple_any(f::Function, t::Tuple) = _tuple_any(f, false, t...)
function _tuple_any(f::Function, tf::Bool, a, b...)
    @_inline_meta
    _tuple_any(f, tf | f(a), b...)
end
_tuple_any(f::Function, tf::Bool) = tf

"""
    empty(x::Tuple)

Returns an empty tuple, `()`.
"""
empty(@nospecialize x::Tuple) = ()

max_stars_repo_name: milesfrain/julia

max_issues_repo_name: AsForMe/julia

max_forks_repo_name: AsForMe/julia


text:
function solve()
    n, a, b = [parse(Int, x) for x in split(readline())]
    n - a + b
end

println(solve())

max_stars_repo_name: aishikawa/atcoder-julia

max_issues_repo_name: aishikawa/atcoder-julia

max_forks_repo_name: aishikawa/atcoder-julia


text:
using SimpleRegressionGP
using Test

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


text:
using DeepDiffs
using Base.Test
using TestSetExtensions

@testset DottedTestSet "DeepDiff Tests" begin
    @includetests ARGS
end


text:
"""
```
Species(id, [value=0])
```

Define a `Species` with a `name` and initial population `value`.
"""
struct Species
  identifier :: Symbol
  population :: Int

  function Species(id, value::Int=0)
    if value < 0
      error("Species population must be nonnegative.")
    end
    new(Symbol(id), value)
  end
end

function Base.show(io::IO, x::Species)
  println(io, x.population)
end


text:
#
# Factors Iterating Tests
#

@testset "Factors Iter" begin

@testset "pattern" begin
    ft = Factor(:l1=>2, :l2=>3)
    @test pattern(ft, :l1) == [1 2 1 2 1 2]'
    @test pattern(ft, [:l1, :l2]) == pattern(ft)
    @test pattern(ft) == [1 1; 2 1; 1 2; 2 2; 1 3; 2 3]
end

@testset "pattern states" begin
    A = Dimension(:A, 2:4)
    B = Dimension(:B, 'a':'e')
    C = Dimension(:C, ["aa", "gg", "ss"])

    ft = Factor([A, B, C])
    len = length(ft)
    manual = Array{Any}(len, 3)
    for (i, a) in enumerate(A)
        for (j, b) in enumerate(B)
            for (k, c) in enumerate(C)
                l = sub2ind(ft, i, j ,k)
                manual[l, :] = [a, b, c]
            end
        end
    end
    @test pattern_states(ft) == manual
end

@testset "pattern iter" begin
    A = Dimension(:A, 2:4)
    B = Dimension(:B, 'a':'e')
    C = Dimension(:C, ["aa", "gg", "ss"])

    ft = Factor([A, B, C])
    len = length(ft)

    for (i, at) in enumerate(ft)
        sub = ind2sub(ft, i)
        @test at == sub2at(ft, sub...)
    end
end

end


text:
using BhAbs_Solver

run_binary_BH_simulation(

	#Time variables
	tspan   =   [0. , 102.]   ,
  	dt                =   0.05        ,
  	out_every         =   1.0         ,
  	max_run_time        =   "24:00:00"  ,

  	#Spatial variables
  	Nnodes              =   128     ,
  	Rcavity            =   40.     ,

  	#Black hole parameters
  	Mbh                =   1.0     ,
  	Rorbit             =   6.0     ,
  	alpha              =   10.0    ,

  	#Initial configuration
  	Gaussian_pulse = [ 3.5 , 3.7 , 20.  , 0.1 , 2 ],
  	)

resume_simulation_from_file(
	"./data/2021-07-27/binary_bh_M_1.0_alpha_10.0_Rorbit_6.0_Rcavity_40.0_N_128_ti_0.0_tf_102.0.h5",
  "./2021-07-27/testfile8.h5",
  t_end = 200.0
  )


text:
#
# Various update methods depending on the solution method
#

# For iterating on initial guess and used within PFI
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{IterateOnPolicy},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid;
    nk, nz = length(kgrid), length(zgrid)

    _inds = LinearIndices((nk, nz))

    # Iterate over all states
    for ik in 1:nk, iz in 1:nz
        # Pull out states
        k = kgrid[ik]
        z = zgrid[iz]

        # Pull out policy and evaluate consumption
        ikiz_index = _inds[ik, iz]
        k1 = kpol[ikiz_index]
        c = expendables_t(ncgm, k, z) - k1

        # New value
        EV = compute_EV(ncgm, vp, k1, iz)
        V[ikiz_index] = u(ncgm, c) + ncgm.β*EV
    end

    # Update coefficients
    update_v!(vp, Φ \ V, 1.0)
    update_k!(vp, Φ \ kpol, 1.0)

    return V
end

# ECM_VFI
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{VFI_ECM},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid;
    nk, nz = length(kgrid), length(zgrid)

    # Iterate over all states
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for ik in 1:nk, iz in 1:nz
        ikiz_index = _inds[ik, iz]
        k = kgrid[ik]
        z = zgrid[iz]

        # Policy from envelope condition
        kp = env_condition_kp!(temp, ncgm, vp, k, z)
        c = expendables_t(ncgm, k, z) - kp
        kpol[ikiz_index] = kp

        # New value
        EV = compute_EV!(temp, ncgm, vp, kp, iz)
        V[ikiz_index] = u(ncgm, c) + ncgm.β*EV
    end

    # Update coefficients
    update_v!(vp, Φ \ V, 1.0)
    update_k!(vp, Φ \ kpol, 1.0)

    return V
end

# VFI_EGM
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{VFI_EGM},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid; grid = ncgm.grid;
    nk, nz = length(kgrid), length(zgrid)

    # Iterate
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for iz=1:nz, ik=1:nk

        # In EGM we use the grid points as if they were our
        # policy for yesterday and find implied kt
        ikiz_index = _inds[ik, iz]
        k1 = kgrid[ik];z = zgrid[iz];

        # Compute the derivative of expected values
        dEV = compute_dEV!(temp, ncgm, vp, k1, iz)

        # Compute optimal consumption
        c = duinv(ncgm, ncgm.β*dEV)

        # Need to find corresponding kt for optimal c
        obj(kt) = expendables_t(ncgm, kt, z) - c - k1
        kt_star = brent(obj, 0.0, 2.0, xtol=1e-10)

        # New value
        EV = compute_EV!(temp, ncgm, vp, k1, iz)
        V[ikiz_index] = u(ncgm, c) + ncgm.β*EV
        kpol[ikiz_index] = kt_star
    end

    # New Φ (has our new "kt_star" and z points)
    Φ_egm = complete_polynomial([kpol grid[:, 2]], vp.d)

    # Update coefficients
    update_v!(vp, Φ_egm \ V, 1.0)
    update_k!(vp, Φ_egm \ grid[:, 1], 1.0)

    # Update V and kpol to be value and policy corresponding
    # to our grid again
    copy!(V, Φ*vp.v_coeffs)
    copy!(kpol, Φ*vp.k_coeffs)

    return V
end

# VFI
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{VFI},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid
    nk, nz = length(kgrid), length(zgrid)

    # Iterate over all states
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for iz=1:nz, ik=1:nk
        k = kgrid[ik]; z = zgrid[iz]

        # Define an objective function (negative for minimization)
        y = expendables_t(ncgm, k, z)
        solme(kp) = du(ncgm, y - kp) - ncgm.β*compute_dEV!(temp, ncgm, vp, kp, iz)

        # Find sol to foc
        kp = brent(solme, 1e-8, y-1e-8; rtol=1e-12)
        c = expendables_t(ncgm, k, z) - kp

        # New value
        ikiz_index = _inds[ik, iz]
        EV = compute_EV!(temp, ncgm, vp, kp, iz)
        V[ikiz_index] = u(ncgm, c) + ncgm.β*EV
        kpol[ikiz_index] = kp
    end

    # Update coefficients
    update_v!(vp, Φ \ V, 1.0)
    update_k!(vp, Φ \ kpol, 1.0)

    return V
end

# PFI ECM
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{PFI_ECM},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Copy valuecoeffs object and use to iterate to
    # convergence given a policy
    vp_igp = copy(vp, IterateOnPolicy())
    solve(ncgm, vp_igp; nskipprint=1000, maxiter=5000, verbose=false)

    # Update the policy and values
    kp = env_condition_kp(ncgm, vp)
    update_k!(vp, Φ \ kp, 1.0)
    update_v!(vp, vp_igp.v_coeffs, 1.0)

    # Update all elements of value
    copy!(V, Φ*vp.v_coeffs)
    copy!(kpol, kp)

    return V
end

# Conventional PFI
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{PFI},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid; grid = ncgm.grid;
    nk, nz = length(kgrid), length(zgrid)

    # Copy valuecoeffs object and use to iterate to
    # convergence given a policy
    vp_igp = copy(vp, IterateOnPolicy())
    solve(ncgm, vp_igp; nskipprint=1000, maxiter=5000, verbose=false)

    # Update the policy and values
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for ik in 1:nk, iz in 1:nz
        k = kgrid[ik]; z = zgrid[iz];

        # Define an objective function (negative for minimization)
        y = expendables_t(ncgm, k, z)
        solme(kp) = du(ncgm, y - kp) - ncgm.β*compute_dEV!(temp, ncgm, vp, kp, iz)

        # Find minimum of objective
        kp = brent(solme, 1e-8, y-1e-8; rtol=1e-12)

        # Update policy function
        ikiz_index = _inds[ik, iz]
        kpol[ikiz_index] = kp
    end

    # Get new coeffs
    update_k!(vp, Φ \ kpol, 1.0)
    update_v!(vp, vp_igp.v_coeffs, 1.0)

    # Update all elements of value
    copy!(V, Φ*vp.v_coeffs)

    return V
end

# Envelope condition iterate on derivative
function update!(dV::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{dVFI_ECM},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    kgrid = ncgm.kgrid; zgrid = ncgm.zgrid; grid = ncgm.grid;
    nk, nz, ns = length(kgrid), length(zgrid), size(grid, 1)

    # Iterate over all states
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for iz=1:nz, ik=1:nk
        k = kgrid[ik]; z = zgrid[iz];

        # Envelope condition implies optimal kp
        kp = env_condition_kp!(temp, ncgm, vp, k, z)
        c = expendables_t(ncgm, k, z) - kp

        # New value
        ikiz_index = _inds[ik, iz]
        dEV = compute_dEV!(temp, ncgm, vp, kp, iz)
        dV[ikiz_index] = (1-ncgm.δ+df(ncgm, k, z))*ncgm.β*dEV
        kpol[ikiz_index] = kp
    end

    # Get new coeffs
    update_k!(vp, Φ \ kpol, 1.0)
    update_v!(vp, dΦ \ dV, 1.0)

    return dV
end

# Conventional Euler equation method
function update!(V::Vector{Float64}, kpol::Vector{Float64},
                 ncgm::NeoclassicalGrowth, vp::ValueCoeffs{EulEq},
                 Φ::Matrix{Float64}, dΦ::Matrix{Float64})
    # Get sizes and allocate for complete_polynomial
    @unpack kgrid, zgrid, weights, z1 = ncgm
    nz1, nz = size(z1)
    nk = length(kgrid)

    # Iterate over all states
    temp = Array{Float64}(undef, n_complete(2, vp.d))
    _inds = LinearIndices((nk, nz))
    for iz in 1:nz, ik in 1:nk
        k = kgrid[ik]; z = zgrid[iz];

        # Create current polynomial
        complete_polynomial!(temp, [k, z], vp.d)

        # Compute what capital will be tomorrow according to policy
        kp = dot(temp, vp.k_coeffs)

        # Compute RHS of EE
        rhs_ee = 0.0
        for iz1 in 1:nz1
            # Possible z in t+1
            zp = z1[iz1, iz]

            # Policy for k_{t+2}
            complete_polynomial!(temp, [kp, zp], vp.d)
            kpp = dot(temp, vp.k_coeffs)

            # Implied t+1 consumption
            cp = expendables_t(ncgm, kp, zp) - kpp

            # Add to running expectation
            rhs_ee += ncgm.β*weights[iz1]*du(ncgm, cp)*(1-ncgm.δ+df(ncgm, kp, zp))
        end

        # The rhs of EE implies consumption and investment in t
        c = duinv(ncgm, rhs_ee)
        kp_star = expendables_t(ncgm, k, z) - c

        # New value
        ikiz_index = _inds[ik, iz]
        EV = compute_EV!(temp, ncgm, vp, kp_star, iz)
        V[ikiz_index] = u(ncgm, c) + ncgm.β*EV
        kpol[ikiz_index] = kp_star
    end

    # Update coefficients
    update_v!(vp, Φ \ V, 1.0)
    update_k!(vp, Φ \ kpol, 1.0)

    return V
end

#
# General solve method
#

function solve(
        ncgm::NeoclassicalGrowth, vp::ValueCoeffs;
        tol::Float64=1e-6, maxiter::Int=5000, dampen::Float64=1.0,
        nskipprint::Int=1, verbose::Bool=true
    )
    # Get number of k and z on grid
    nk, nz = length(ncgm.kgrid), length(ncgm.zgrid)

    # Build basis matrix and value function
    dPhi = complete_polynomial(ncgm.grid, vp.d, Derivative{1}())
    Phi = complete_polynomial(ncgm.grid, vp.d)
    V = build_V_or_dV(ncgm, vp)
    k = build_k(ncgm, vp)
    Vnew = copy(V)
    knew = copy(k)

    # Print column names
    if verbose
        @printf("| Iteration | Distance V | Distance K |\n")
    end

    # Iterate to convergence
    dist, iter = 10.0, 0
    while (tol < dist) & (iter < maxiter)
        # Update the value function using appropriate update methods
        update!(Vnew, knew, ncgm, vp, Phi, dPhi)

        # Compute distance and update all relevant elements
        iter += 1
        dist_v = maximum(abs, 1.0 .- Vnew./V)
        dist_k = maximum(abs, 1.0 .- knew./k)
        copy!(V, Vnew)
        copy!(k, knew)

        # If we are iterating on a policy, use the difference of values
        # otherwise use the distance on policy
        dist = ifelse(solutionmethod(vp) == IterateOnPolicy, dist_v, dist_k)

        # Print status update
        if verbose && (iter%nskipprint == 0)
            @printf("|%-11d|%-12e|%-12e|\n", iter, dist_v, dist_k)
        end
    end

    # Update value and policy functions one last time as long as the
    # solution method isn't IterateOnPolicy
    if ~(solutionmethod(vp) == IterateOnPolicy)
        # Update capital policy after finished
        kp = env_condition_kp(ncgm, vp)
        update_k!(vp, complete_polynomial(ncgm.grid, vp.d) \ kp, 1.0)

        # Update value function according to specified policy
        vp_igp = copy(vp, IterateOnPolicy())
        solve(ncgm, vp_igp; tol=1e-10, maxiter=5000, verbose=false)
        update_v!(vp, vp_igp.v_coeffs, 1.0)

    end

    return vp
end

max_stars_repo_name: sglyon/CLMMJuliaPythonMatlab

max_issues_repo_name: sglyon/CLMMJuliaPythonMatlab

max_forks_repo_name: sglyon/CLMMJuliaPythonMatlab


text:
module RemoteSemaphores

export RemoteSemaphore, acquire, release

using Base: Semaphore, acquire, release
using Distributed

"""
    RemoteSemaphore(n::Int, pid=myid())

A semaphore living on a specific process.
Do not attempt to fetch the future to a different process and use it there, as that will be
an isolated, unsynced copy of the semaphore.
"""
struct RemoteSemaphore
    n::Int  # stored for printing only
    rref::Future

    function RemoteSemaphore(n::Integer, pid=myid())
        sem = Semaphore(n)
        fut = Future(pid)
        put!(fut, sem)

        return new(n, fut)
    end
end

function Base.acquire(rsem::RemoteSemaphore)
    fut = rsem.rref
    loc = fut.where
    remotecall_wait(loc) do
        sem = fetch(fut)
        acquire(sem)
    end

    return nothing
end

function Base.release(rsem::RemoteSemaphore)
    fut = rsem.rref
    loc = fut.where
    remotecall_wait(loc) do
        sem = fetch(fut)
        release(sem)
    end

    return nothing
end

function Base.show(io::IO, rsem::RemoteSemaphore)
    print(io, typeof(rsem), '(', rsem.n, ", pid=", rsem.rref.where, ')')
end

# expensive, easily out of sync, for testing only
function _current_count(rsem::RemoteSemaphore)
    fut = rsem.rref
    loc = fut.where

    return remotecall_fetch(loc) do
        sem = fetch(fut)
        return sem.curr_cnt
    end
end

end


text:
"""
    BasisFunctionRegressor{Tblr,Tϕ}

A Basis Function Regressor represents a Bayesian Linear Regressor where the
input `x` is first mapped to a feature space through a basis function ϕ.

ϕ must be a function which accepts one of the allowed input types for
BayesianLinearRegressors (ColVecs, RowVecs or Matrix{<:Real} - see the package
readme for more details) and it must output one of these allowed types.

```jldoctest; output = false
x = RowVecs(hcat(range(-1.0, 1.0, length=5)))
blr = BayesianLinearRegressor(zeros(2), Diagonal(ones(2)))

ϕ(x::RowVecs) = RowVecs(hcat(ones(length(x)), prod.(x)))
bfr = BasisFunctionRegressor(blr, ϕ)

var(bfr(x))

# output

5-element Vector{Float64}:
 2.0
 1.25
 1.0
 1.25
 2.0
```

See [1], Section 3.1 for more details on basis function regression.

[1] - C. M. Bishop. "Pattern Recognition and Machine Learning". Springer, 2006.
"""
struct BasisFunctionRegressor{Tblr<:BayesianLinearRegressor,Tϕ} <: AbstractGP
    blr::Tblr
    ϕ::Tϕ
end

function (bfblr::BasisFunctionRegressor)(x::AbstractVector, args...)
    return bfblr.blr(bfblr.ϕ(x), args...)
end


text:
module Tilings

using Luxor
using Colors

include("tiles.jl")

export PolyTile, draw

end


text:
function applyliverics!(system::CVSystem,old=Dict("a"=>0),restart="no")
    if restart == "no"
        system.liver.P[1,1] = 4.99*mmHgToPa;
        system.liver.V[1,1] = (system.liver.P[1,1]*system.liver.C[1] +
            system.liver.V0[1]);
        system.liver.P[1,2] = 4.91*mmHgToPa;
        system.liver.V[1,2] = (system.liver.P[1,2]*system.liver.C[2] +
            system.liver.V0[2]);
        system.liver.P[1,3] = 4.79*mmHgToPa;
        system.liver.V[1,3] = (system.liver.P[1,3]*system.liver.C[3] +
            system.liver.V0[3]);
        system.liver.P[1,4] = 4.66*mmHgToPa;
        system.liver.V[1,4] = (system.liver.P[1,4]*system.liver.C[4] +
            system.liver.V0[4]);
    elseif restart == "yes"
        system.liver.P[1,1] = old["P"][end,1];
        system.liver.V[1,1] = old["V"][end,1];
        system.liver.Q[1,1] = old["Q"][end,1];
        system.liver.P[1,2] = old["P"][end,2];
        system.liver.V[1,2] = old["V"][end,2];
        system.liver.Q[1,2] = old["Q"][end,2];
        system.liver.P[1,3] = old["P"][end,3];
        system.liver.V[1,3] = old["V"][end,3];
        system.liver.Q[1,3] = old["Q"][end,3];
        system.liver.P[1,4] = old["P"][end,4];
        system.liver.V[1,4] = old["V"][end,4];
        system.liver.Q[1,4] = old["Q"][end,4];
    end
end

max_stars_repo_name: dcanuto/CardioModeling

max_issues_repo_name: dcanuto/CardioModeling

max_forks_repo_name: dcanuto/CardioModeling


text:
using DataFrames
dvector = data([1,2,3,4,5])
dmatrix = data([10 20 30; 40 50 60])

max_stars_repo_path: Chapter 2/5.jl

max_stars_repo_name: PacktPublishing/Julia-for-data-science

max_issues_repo_path: Chapter 2/5.jl

max_issues_repo_name: PacktPublishing/Julia-for-data-science

max_forks_repo_path: Chapter 2/5.jl

max_forks_repo_name: PacktPublishing/Julia-for-data-science


text:
# Surrogate Model with Gaussian Processes
module GPs

using ..DE

using LinearAlgebra

export model_training, predict

function kernel_trick(x1::Float64, x2::Float64,
                      σ_f::Float64, l::Float64)
    re = σ_f^2 * exp((x1 - x2)^2 / -2)
end

function K(x::Array{Float64}, σ_f::Float64, l::Float64)
    le = length(x)
    re = Array{Float64}(undef, le, le)
    for i=1:1:le
        for j=1:1:le
            re[i, j] = kernel_trick(x[i], x[j], σ_f, l)
        end
    end
    re
end

function K(x1::Array{Float64}, x2::Float64, σ_f::Float64, l::Float64)
    le = length(x1)
    re = Array{Float64}(undef, le, 1)
    for i=1:1:le
        re[i, 1] = kernel_trick(x1[i], x2, σ_f, l)
    end
    re
end

function K(x::Float64, σ_f::Float64, l::Float64)
    kernel_trick(x, x, σ_f, l)
end

function marginal_likelihood(tr_x::Array{Float64}, tr_y::Array{Float64},
                             σ_f::Float64, l::Float64, σ_n::Float64)
    A = K(tr_x, σ_f, l) + σ_n^2 * UniformScaling(length(tr_x))
    n = length(tr_x)
    -(-1/2 * tr_y' * inv(A) * tr_y - 1/2 * log(det(A)) - n/2 * log(2 * pi))
end

function model_training(tr_x::Array{Float64}, tr_y::Array{Float64})
    f(pop::DE.Population) = begin
        for i in 1:size(pop.x, 1)
            pop.fit[i] = marginal_likelihood(tr_x, tr_y,
                            pop.x[i, 1], pop.x[i, 2], pop.x[i, 3])
        end
    end

    σ_f, l, σ_n= DE_evolution(f, 10, 3, 50)[1]
    A = K(tr_x, σ_f, l) + σ_n * UniformScaling(length(tr_x))
    while !isposdef(A)
        σ_f, l, σ_n= DE_evolution(f, 10, 3, 50)[1]
        A = K(tr_x, σ_f, l) + σ_n * UniformScaling(length(tr_x))
    end

    L = cholesky(A)
    α = L.U\(L.L\tr_y)
    return σ_f, l, L, α
end

function predict(σ_f::Float64, l::Float64, L::Cholesky, α::Array{Float64},
                 tr_x::Array{Float64}, pr_x::Float64)
    k_s = K(tr_x, pr_x, σ_f, l)
    k_ss = K(pr_x, σ_f, l)
    mean = k_s' ⋅ α
    v = L\k_s
    vari = k_ss - v' ⋅ v
    return mean, vari
end

end  # modue SuM


text:
## wrapper for librdkafka C API, see:
## https://github.com/edenhill/librdkafka/blob/master/src/rdkafka.h

function rd_kafka_version()
    ccall((:rd_kafka_version, RDKafka.librdkafka), Cint, ())
end

## rd_kafka_conf_t

function kafka_conf_new()
    return ccall((:rd_kafka_conf_new, librdkafka), Ptr{Cvoid}, ())
end

function kafka_conf_destroy(conf::Ptr{Cvoid})
    ccall((:rd_kafka_conf_destroy, librdkafka), Cvoid, (Ptr{Cvoid},), conf)
end

function kafka_conf_set(conf::Ptr{Cvoid}, key::String, val::String)
    err_str = zeros(UInt8, 512)
    return ccall((:rd_kafka_conf_set, librdkafka), Cvoid,
                 (Ptr{Cvoid}, Cstring, Cstring, Ptr{UInt8}, Csize_t),
                 conf, key, val, pointer(err_str), sizeof(err_str))
end

function kafka_conf_get(conf::Ptr{Cvoid}, key::String)
    dest = zeros(UInt8, 512)
    sz_ref = Ref{Csize_t}(0)
    ccall((:rd_kafka_conf_get, librdkafka), Cvoid,
          (Ptr{Cvoid}, Cstring, Ptr{UInt8}, Ptr{Csize_t}),
          conf, key, pointer(dest), sz_ref)
    return unsafe_string(pointer(dest))
end

function kafka_conf_set_error_cb(conf::Ptr{Cvoid}, c_fn::Ptr{Cvoid})
    ccall((:rd_kafka_conf_set_error_cb, librdkafka), Cvoid,
          (Ptr{Cvoid}, Ptr{Cvoid}), conf, c_fn)
end

function kafka_conf_set_dr_msg_cb(conf::Ptr{Cvoid}, c_fn::Ptr{Cvoid})
    ccall((:rd_kafka_conf_set_dr_msg_cb, librdkafka), Cvoid,
          (Ptr{Cvoid}, Ptr{Cvoid}), conf, c_fn)
end

## rd_kafka_t

const KAFKA_TYPE_PRODUCER = Cint(0)
const KAFKA_TYPE_CONSUMER = Cint(1)

function kafka_new(conf::Ptr{Cvoid}, kafka_type::Cint)
    err_str = zeros(UInt8, 512)
    client = ccall((:rd_kafka_new, librdkafka),
                   Ptr{Cvoid},
                   (Cint, Ptr{Cvoid}, Ptr{UInt8}, Csize_t),
                   kafka_type, conf, pointer(err_str), sizeof(err_str))
    return client
end

function kafka_destroy(rk::Ptr{Cvoid})
    ccall((:rd_kafka_destroy, librdkafka), Cvoid, (Ptr{Cvoid},), rk)
end

function kafka_consumer_close(rk::Ptr{Cvoid})
    ccall((:rd_kafka_consumer_close, librdkafka), Cvoid, (Ptr{Cvoid},), rk)
end

## rd_kafka_topic_conf_t

function kafka_topic_conf_new()
    return ccall((:rd_kafka_topic_conf_new, librdkafka), Ptr{Cvoid}, ())
end

function kafka_topic_conf_destroy(conf::Ptr{Cvoid})
    ccall((:rd_kafka_topic_conf_destroy, librdkafka), Cvoid, (Ptr{Cvoid},), conf)
end

function kafka_topic_conf_set(conf::Ptr{Cvoid}, key::String, val::String)
    err_str = Array{UInt8}(undef, 512)
    return ccall((:rd_kafka_topic_conf_set, librdkafka), Cvoid,
                 (Ptr{Cvoid}, Cstring, Cstring, Ptr{UInt8}, Csize_t),
                 conf, key, val, pointer(err_str), sizeof(err_str))
end

function kafka_topic_conf_get(conf::Ptr{Cvoid}, key::String)
    dest = Array{UInt8}(undef, 512)
    sz_ref = Ref{Csize_t}(0)
    ccall((:rd_kafka_topic_conf_get, librdkafka), Cvoid,
          (Ptr{Cvoid}, Cstring, Ptr{UInt8}, Ptr{Csize_t}),
          conf, key, pointer(dest), sz_ref)
    return unsafe_string(pointer(dest))
end

# rd_kafka_topic_t

function kafka_topic_new(rk::Ptr{Cvoid}, topic::String, topic_conf::Ptr{Cvoid})
    return ccall((:rd_kafka_topic_new, librdkafka), Ptr{Cvoid},
                 (Ptr{Cvoid}, Cstring, Ptr{Cvoid}),
                 rk, topic, topic_conf)
end

function kafka_topic_destroy(rkt::Ptr{Cvoid})
    ccall((:rd_kafka_topic_destroy, librdkafka), Cvoid, (Ptr{Cvoid},), rkt)
end

function kafka_poll(rk::Ptr{Cvoid}, timeout::Integer)
    return ccall((:rd_kafka_poll, librdkafka), Cint,
                 (Ptr{Cvoid}, Cint),
                 rk, timeout)

end

function produce(rkt::Ptr{Cvoid}, partition::Integer,
                 key::Vector{UInt8}, payload::Vector{UInt8})
    flags = Cint(0)
    errcode = ccall((:rd_kafka_produce, librdkafka), Cint,
                    (Ptr{Cvoid}, Int32, Cint,
                     Ptr{Cvoid}, Csize_t,
                     Ptr{Cvoid}, Csize_t,
                     Ptr{Cvoid}),
                    rkt, Int32(partition), flags,
                    pointer(payload), length(payload),
                    pointer(key), length(key),
                    C_NULL)

    if errcode != 0
        #error("Produce request failed with error code (unix): $errcode")
        ## errno is c global var
	errnr = unsafe_load(cglobal(:errno, Int32))
        errmsg = unsafe_string(ccall(:strerror, Cstring, (Int32,), errnr))
        error("Produce request failed with error code: $errmsg")
    end
end

## topic partition list

mutable struct CKafkaTopicPartition
    topic::Cstring
    partition::Cint
    offset::Clonglong
    metadata::Ptr{Cvoid}
    metadata_size::Csize_t
    opaque::Ptr{Cvoid}
    err::Cint
    _private::Ptr{Cvoid}
end
mutable struct CKafkaTopicPartitionList
    cnt::Cint
    size::Cint
    elems::Ptr{CKafkaTopicPartition}
end

function kafka_topic_partition_list_new(sz::Integer=0)
    rkparlist = ccall((:rd_kafka_topic_partition_list_new, librdkafka), Ptr{CKafkaTopicPartitionList},
                       (Cint,), sz)
    if rkparlist != Ptr{CKafkaTopicPartitionList}(0)
        return rkparlist
    end
    return nothing
end

function kafka_topic_partition_list_destroy(rkparlist::Ptr{CKafkaTopicPartitionList})
    ccall((:rd_kafka_topic_partition_list_destroy, librdkafka), Cvoid, (Ptr{CKafkaTopicPartitionList},), rkparlist)
end

function kafka_topic_partition_list_add(rkparlist::Ptr{CKafkaTopicPartitionList},
                                        topic::String, partition::Integer)
    ccall((:rd_kafka_topic_partition_list_add, librdkafka), Ptr{CKafkaTopicPartition},
          (Ptr{CKafkaTopicPartitionList}, Cstring, Int32,), rkparlist, topic, partition)
end

## partition assignment

function kafka_assignment(rk::Ptr{Cvoid}, rkparlist::Ptr{CKafkaTopicPartitionList})
    errcode = ccall((:rd_kafka_assignment, librdkafka), Cint,
                    (Ptr{Cvoid}, Ref{Ptr{CKafkaTopicPartitionList}}), rk, rkparlist)
    if errcode != 0
        error("Assignment retrieval failed with error $errcode")
    end
end

## subscribe

function kafka_subscribe(rk::Ptr{Cvoid}, rkparlist::Ptr{CKafkaTopicPartitionList})
    errcode = ccall((:rd_kafka_subscribe, librdkafka), Cint,
                    (Ptr{Cvoid}, Ptr{Cvoid}), rk, rkparlist)
    if errcode != 0
        error("Subscription failed with error $errcode")
    end
end

## seek

function kafka_offsets_for_times(rk::Ptr{Cvoid}, rkparlist::Ptr{CKafkaTopicPartitionList}, timeout::Integer)
    errcode = ccall((:rd_kafka_offsets_for_times, librdkafka), Cint,
                    (Ptr{Cvoid}, Ptr{CKafkaTopicPartitionList}, Cint), rk, rkparlist, timeout)
    if errcode != 0
        error("kafka_offsets_for_times failed with error $errcode")
    end
end

function kafka_seek(rkt::Ptr{Cvoid}, partition::Int32, offset::Int64, timeout::Int=1000)
    errcode = ccall((:rd_kafka_seek, librdkafka), Cint,
                    (Ptr{Cvoid}, Cint, Clonglong, Cint),
                    rkt, partition, offset, timeout)
    if errcode != 0
        error("kafka_seek failed with error $errcode")
    end
end

function kafka_assign(rk::Ptr{Cvoid}, rkparlist::Ptr{CKafkaTopicPartitionList})
    errcode = ccall((:rd_kafka_assign, librdkafka), Cint,
                    (Ptr{Cvoid}, Ptr{CKafkaTopicPartitionList}),
                    rk, rkparlist)
    if errcode != 0
        error("kafka_assign failed with error $errcode")
    end
end

struct CKafkaMessage
    err::Cint
    rkt::Ptr{Cvoid}
    partition::Int32
    payload::Ptr{UInt8}
    len::Csize_t
    key::Ptr{UInt8}
    key_len::Csize_t
    offset::Int64
    _private::Ptr{Cvoid}
end

function kafka_consumer_poll(rk::Ptr{Cvoid}, timeout::Integer)
    msg_ptr = ccall((:rd_kafka_consumer_poll, librdkafka), Ptr{CKafkaMessage},
                    (Ptr{Cvoid}, Cint), rk, timeout)
    if msg_ptr != Ptr{CKafkaMessage}(0)
        return msg_ptr
    else
        return nothing
    end

end

function kafka_message_destroy(msg_ptr::Ptr{CKafkaMessage})
    ccall((:rd_kafka_message_destroy, librdkafka), Cvoid, (Ptr{Cvoid},), msg_ptr)
end


text:
#==============================================================================#
# MobileAnalytics.jl
#
# This file is generated from:
# https://github.com/aws/aws-sdk-js/blob/master/apis/mobileanalytics-2014-06-05.normal.json
#==============================================================================#

__precompile__()

module MobileAnalytics

using AWSCore

"""
    using AWSSDK.MobileAnalytics.put_events
    put_events([::AWSConfig], arguments::Dict)
    put_events([::AWSConfig]; events=, *header:* x-amz-Client-Context=, <keyword arguments>)

    using AWSCore.Services.mobileanalytics
    mobileanalytics([::AWSConfig], "POST", "/2014-06-05/events", arguments::Dict)
    mobileanalytics([::AWSConfig], "POST", "/2014-06-05/events", events=, *header:* x-amz-Client-Context=, <keyword arguments>)

# PutEvents Operation

The PutEvents operation records one or more events. You can have up to 1,500 unique custom events per app, any combination of up to 40 attributes and metrics per custom event, and any number of attribute or metric values.

# Arguments

## `events = [[ ... ], ...]` -- *Required*
An array of Event JSON objects
```
 events = [[
        "eventType" => <required> ::String,
        "timestamp" => <required> ::String,
        "session" =>  [
            "id" =>  ::String,
            "duration" =>  ::Int,
            "startTimestamp" =>  ::String,
            "stopTimestamp" =>  ::String
        ],
        "version" =>  ::String,
        "attributes" =>  ::Dict{String,String},
        "metrics" =>  ::Dict{String,String}
    ], ...]
```

## `*header:* x-amz-Client-Context = ::String` -- *Required*
The client context including the client ID, app title, app version and package name.

## `*header:* x-amz-Client-Context-Encoding = ::String`
The encoding used for the client context.

# Exceptions

`BadRequestException`.

See also: [AWS API Documentation](https://docs.aws.amazon.com/goto/WebAPI/mobileanalytics-2014-06-05/PutEvents)
"""
@inline put_events(aws::AWSConfig=default_aws_config(); args...) = put_events(aws, args)

@inline put_events(aws::AWSConfig, args) = AWSCore.Services.mobileanalytics(aws, "POST", "/2014-06-05/events", args)

@inline put_events(args) = put_events(default_aws_config(), args)

end # module MobileAnalytics

#==============================================================================#
# End of file
#==============================================================================#