text:
#
# Calculates requests per second from client to server without response
# Results: rps: 43796.25990550989, time: 22.83299994468689
#
include("../../src/global/ImportFolders.jl")

import Config
import Server
import Connection
import Event

function onBeforeResponse(sock::Base.TCPSocket, data::Array{Any, 1}, ans::Connection.Answer)
  for i=1:1000000 Server.request(sock, UInt8(1), UInt(123)) end
end

con = Server.create(ip"127.0.0.1", Config.val(:conServerPort), true)
Event.on(con.observer, Connection.EVENT_BEFORE_RESPONSE, onBeforeResponse)
Server.run(con)

max_stars_repo_name: tmptrash/jevo

max_issues_repo_name: tmptrash/jevo

max_forks_repo_name: tmptrash/jevo


text:
using IntervalUnionArithmetic
using Test

include("interval_unions.jl")
include("set_operations.jl")
include("arithmetic.jl")

ext: jl


text:
using Serialization

macro gctime(ex)
    fc = isdefined(Base.Experimental, Symbol("@force_compile")) ?
        :(Base.Experimental.@force_compile) :
        :()
    quote
        $fc
        local start_gc_num = Base.gc_num()
        local start_time = time_ns()
        local val = $(esc(ex))
        local end_time = time_ns()
        local end_gc_num = Base.gc_num()
        result = (
            value = val,
            times = (end_time - start_time),
            stats = Base.GC_Diff(end_gc_num, start_gc_num),
        )
        "SERIALIZE" in ARGS ? serialize(stdout, result) : display(result)
    end
end

max_stars_repo_name: oscardssmith/GCBenchmarks

max_issues_repo_name: oscardssmith/GCBenchmarks

max_forks_repo_name: oscardssmith/GCBenchmarks


text:
# The goal of these methods is to determine the JCF of the matrix
# of tree functions.

using DiscreteFunctions, LinearAlgebra, Counters, SimpleTools

function find_nonzero(v::Vector{T}) where T<:Real
    idx = findall(abs.(v) .> 0.1)
    return idx[1]
end

function cols2nums(NS::Matrix{T}) where T<:Real
    r,c = size(NS)
    idxlist = [ find_nonzero(NS[:,j]) for j=1:c ]
    return Set{Int}(idxlist)
end

"""
`null_vecs(A)`: We assume the nullspace of `A` is spanned by standard
basis vectors. This function determines which standard basis vectors
those are. The return is a set of integers.
"""
function null_vecs(A::Matrix{T})::Set{Int} where T<:Real
    NS = nullspace(A)
    return cols2nums(NS)
end

"""
`generalized_null(A)` computes the nullspace of `A`, `A^2`, `A^3`, and
so on. It returns the standard basis vectors of those nullspaces. The
first set is the nullspace of `A`. The second set contains the additional
vectors found in the nullspace of `A^2`, and so on.

`generalized_null(f)` is a short cut for `generalized_null(Matrix(f))`
"""
function generalized_null(A::Matrix{T}) where T<:Real
    found = Set{Int}()
    list  = Set{Int}[]
    k=1
    while true
        vecs = null_vecs(A^k)
        if length(vecs) == length(found)  # nothing new, we're done
            return list
        end
        newvecs = vecs
        push!(list, vecs)
        found = vecs
        k += 1
    end
end

generalized_null(f::DiscreteFunction) = generalized_null(Matrix(f))

"""
`jordan(A::Matrix)` or `jordan(f::DiscreteFunction)` reports the Jordan structure
for the zero eigenvalues of `A` or `Matrix(f)`.
"""
function jordan(A::Matrix{T}) where T<:Real
    sets = generalized_null(A)
    n = length(sets)
    r = length.(sets)
    s = 0r
    s[1] = r[1]
    for j=2:n
        s[j] = r[j]-r[j-1]
    end
    m = 0s

    for j=1:n-1
        m[j] = s[j]-s[j+1]
    end
    m[n] = s[n]

    blocks = Counter{Int}()
    for j=1:n
        blocks[j] += m[j]
    end

    return blocks
end

jordan(f::DiscreteFunction) = jordan(Matrix(f))

function _zero_block(k::Int)
    A = zeros(Int,k,k)
    for j=1:k-1
        A[j,j+1]=1
    end
    return A
end

"""
`jordan_report(f::DiscreteFunction)` provides information about
the JCF of `Matrix(f)` where we assume `f` is a tree function.
The JCF is returned. We assume `f` has exactly one eigenvalue equal
to 1 and the rest of the eigenvalues are 0s.
"""
function jordan_report(A::Matrix{T}) where T<:Real
    counts = jordan(A)
    nc = length(counts)
    for j=1:nc
        if counts[j]>0
            println("$(counts[j]) Jordan blocks of size $j")
        end
    end

    JCF = ones(Int,1,1)
    for j=nc:-1:1
        for a = 1:counts[j]
            JCF = dcat(JCF,_zero_block(j))
        end
    end

    return JCF
end
jordan_report(f::DiscreteFunction) = jordan_report(Matrix(f))


text:
# --------------
# Initialization
# --------------
using MUMIMO

# --------------
# Basic test
# --------------
mimoUplink(100,4,[4],[4],[0],[12],12,[-5.0:5:20;])

max_stars_repo_name: christianpeel/MUMIMO

max_issues_repo_name: christianpeel/MUMIMO

max_forks_repo_name: christianpeel/MUMIMO


text:
using Test
import Bresenham

function test_line(args...)
  points = Vector{Tuple{Int, Int}}()
  Bresenham.line(args...) do x, y
    push!(points, (x, y))
  end
  points
end

# Low slope (X major)
@test test_line(1, 1, 5, 3) == [ (1, 1), (2, 1), (3, 2), (4, 2), (5, 3) ]
# Reverse low slope (X major)
@test test_line(5, 3, 1, 1) == reverse!([ (1, 1), (2, 1), (3, 2), (4, 2), (5, 3) ])

# Oblique slope (X major)
@test test_line(1, 1, 5, 5) == [ (1, 1), (2, 2), (3, 3), (4, 4), (5, 5) ]
# Reverse oblique slope (X major)
@test test_line(5, 5, 1, 1) == reverse!([ (1, 1), (2, 2), (3, 3), (4, 4), (5, 5) ])

# High slope (Y major)
@test test_line(1, 1, 3, 5) == [ (1, 1), (1, 2), (2, 3), (2, 4), (3, 5) ]
# Reverse high slope (Y major)
@test test_line(3, 5, 1, 1) == reverse!([ (1, 1), (1, 2), (2, 3), (2, 4), (3, 5) ])

# Horizontal
@test test_line(1, 1, 5, 1) == [ (1, 1), (2, 1), (3, 1), (4, 1), (5, 1) ]
# Reverse horizontal
@test test_line(5, 1, 1, 1) == reverse!([ (1, 1), (2, 1), (3, 1), (4, 1), (5, 1) ])
# Vertical
@test test_line(1, 1, 1, 5) == [ (1, 1), (1, 2), (1, 3), (1, 4), (1, 5) ]
# Reverse vertical
@test test_line(1, 5, 1, 1) == reverse!([ (1, 1), (1, 2), (1, 3), (1, 4), (1, 5) ])


text:
"""
update_velocities!(model::AbstractModel)

Solve momentum equation to update the velocities, plus Picard iteration for non-linear rheology.

"""
function update_velocities!(model::AbstractModel)
@unpack params,solver_params=model
@unpack gu,gv,wu,wv = model.fields
n = gu.n + gv.n

x=get_start_guess(model)
b=get_rhs(model)

rel_resid=zero(eltype(b))

converged::Bool = false
i_picard::Int64 = 0
while !converged && (i_picard < solver_params.maxiter_picard)

    i_picard = i_picard + 1

    set_velocities!(model,x)
    update_shelf_strain_rate!(model)
    update_av_speed!(model)
    update_bed_speed!(model)
    update_β!(model)
    update_basal_drag!(model)
    inner_update_viscosity!(model)
    update_av_viscosity!(model)
    update_quadrature_falpha!(model)
    update_βeff!(model)
    update_βeff_on_uv_grids!(model)
    update_rheological_operators!(model)
    op=get_op(model)

    rel_resid = norm(b .- op*x)/norm(b)
    converged = rel_resid < solver_params.tol_picard

    p=get_preconditioner(model,op)
    precondition!(x, p, b)

end
set_velocities!(model,x)

return model
end

"""
    get_start_guess(model::AbstractModel)

 Return starting guess used to begin iterative solution of velocities.

"""
function get_start_guess(model::AbstractModel)
    @unpack gu,gv=model.fields
    @assert eltype(gu.u)==eltype(gv.v)
    n = gu.n + gv.n
    x=[gu.samp*gu.u[:];gv.samp*gv.v[:]]
    return x
end

"""
    get_rhs(model::AbstractModel)

 Return right hand side vector of momentum equations.

"""
function get_rhs(model::AbstractModel)
    @unpack gh,gu,gv,gc=model.fields
    @unpack params=model
    onesvec=ones(gh.nxh*gh.nyh)
    surf_elev_adjusted = gh.crop*(gh.s[:] .+ params.dt*gh.dsdh[:].*(gh.accumulation[:].-gh.basal_melt[:]))
    f1=[
        (params.density_ice*params.g*gu.h[gu.mask]).*(gu.samp*(-gu.∂x'*surf_elev_adjusted))
        ;
        (params.density_ice*params.g*gv.h[gv.mask]).*(gv.samp*(-gv.∂y'*surf_elev_adjusted))
       ]
    f2=[
        (0.5*params.density_ice*params.g*gu.h[gu.mask].^2
        -0.5*params.density_ocean*params.g*(icedraft.(gu.s[gu.mask],gu.h[gu.mask],params.sea_level_wrt_geoid)).^2
        -params.density_ice*params.g*gu.h[gu.mask].*gu.s[gu.mask]).*gu.samp*(-gu.∂x'*(gh.crop*onesvec))
        ;
        (0.5*params.density_ice*params.g*gv.h[gv.mask].^2
        -0.5*params.density_ocean*params.g*(icedraft.(gv.s[gv.mask],gv.h[gv.mask],params.sea_level_wrt_geoid)).^2
        -params.density_ice*params.g*gv.h[gv.mask].*gv.s[gv.mask]).*gv.samp*(-gv.∂y'*(gh.crop*onesvec))
        ]
    rhs=f1+f2
    return rhs
end

"""
    set_velocities!(model::AbstractModel,x)

Set velocities to particular values. Input vector x represents stacked u and v components at valid grid points.
"""
function set_velocities!(model::AbstractModel,x)
    @unpack gh,gu,gv,gc=model.fields
    gu.u[:] .= gu.spread*x[1:gu.n]
    gv.v[:] .= gv.spread*x[(gu.n+1):(gu.n+gv.n)]
    return model
end

"""
    update_shelf_strain_rate!(model::AbstractModel)

Find the effective strain rate for 'ice shelf' parts of strain rate tensor, neglecting all vertical shear.
"""
function update_shelf_strain_rate!(model::AbstractModel)
    @unpack gh,gu,gv,gc=model.fields
    gh.shelf_strain_rate[:] .= sqrt.( (gh.crop*(gu.∂x*(gu.crop*gu.u[:]))).^2 .+
                                      (gh.crop*(gv.∂y*(gv.crop*gv.v[:]))).^2 .+
                                (gh.crop*(gu.∂x*(gu.crop*gu.u[:]))).*(gh.crop*(gv.∂y*(gv.crop*gv.v[:]))) .+
                       0.25*(gh.crop*(gc.cent*(gc.crop*( gu.∂y*(gu.crop*gu.u[:]) .+ gv.∂x*(gv.crop*gv.v[:]) )))).^2  )
    return model
end

"""
    update_av_speed!(model::AbstractModel)

Find the depth-averaged speed on the h-grid using components on u- and v- grids
"""
function update_av_speed!(model::AbstractModel)
    @unpack gh,gu,gv=model.fields
    gh.av_speed[:] .= sqrt.( (gh.crop*(gu.cent*(gu.crop*gu.u[:]))).^2 .+ (gh.crop*(gv.cent*(gv.crop*gv.v[:]))).^2 )
    return model
end

"""
    update_bed_speed!(model::AbstractModel)

Find the sliding speed at the bed on the h-grid using the average speed.
"""
function update_bed_speed!(model::AbstractModel)
    @unpack gh=model.fields
    gh.bed_speed .= gh.av_speed ./ (1.0 .+ gh.quad_f2 .* gh.β)
    return model
end

"""
    update_β!(model::AbstractModel)

Find the drag coefficient at the bed using the sliding law.
"""
function update_β!(model::AbstractModel)
    @unpack gh=model.fields
    @unpack params=model
    gh.β .= gh.weertman_c .* ( sqrt.(gh.bed_speed.^2 .+  params.weertman_reg_speed^2 ) ).^(1.0/params.weertman_m - 1.0)
    return model
end

"""
    update_basal_drag!(model::AbstractModel)

Find the shear stress at the bed.
"""
function update_basal_drag!(model::AbstractModel)
    @unpack gh=model.fields
    gh.τbed .= gh.β .* gh.bed_speed
    return model
end

"""
    inner_update_viscosity!(model::AbstractModel)

Inner update to iteratively refine viscosity on the 3d grid at all sigma levels.
"""
function inner_update_viscosity!(model::AbstractModel)
    @unpack gh,g3d=model.fields
    @unpack params,solver_params=model
    for k=1:g3d.nσs
        for j=1:g3d.nys
            for i=1:g3d.nxs
                if gh.mask[i,j]
                    for iter=1:solver_params.n_iter_viscosity
                        g3d.η[i,j,k] = 0.5 * g3d.glen_b[i,j,k] * (
                                                   sqrt(    gh.shelf_strain_rate[i,j]^2 +
                                                            0.25*(gh.τbed[i,j]*g3d.ζ[k]/g3d.η[i,j,k])^2 +
                                                            params.glen_reg_strain_rate^2   )
                                                                 )^(1.0/params.glen_n - 1.0)
                    end
                end
            end
        end
    end
    return model
end

"""
    update_av_viscosity!(model::AbstractModel)

Use quadrature to compute the depth averaged viscosity.
"""
function update_av_viscosity!(model::AbstractModel)
    @unpack gh,g3d=model.fields
    gh.ηav .= zero(gh.ηav)
    for k=1:g3d.nσs
       for j = 1:g3d.nys
          for i = 1:g3d.nxs
                gh.ηav[i,j] += g3d.quadrature_weights[k] * g3d.η[i,j,k]
          end
       end
    end
    return model
end

"""
    update_quadrature_falpha!(model::AbstractModel)

Use quadrature to compute falpha functions, used to relate average velocities, basal velocities, and surface velocities to one another
"""
function update_quadrature_falpha!(model::AbstractModel)
    @unpack gh,g3d=model.fields
    gh.quad_f1 .= zero(gh.quad_f1)
    gh.quad_f2 .= zero(gh.quad_f2)
    for k=1:g3d.nσs
       for j = 1:g3d.nys
          for i = 1:g3d.nxs
                gh.quad_f1[i,j] += g3d.quadrature_weights[k]*gh.h[i,j]*g3d.ζ[k]/g3d.η[i,j,k]
                gh.quad_f2[i,j] += g3d.quadrature_weights[k]*gh.h[i,j]*(g3d.ζ[k])^2/g3d.η[i,j,k]
          end
       end
    end
    return model
end

"""
    update_βeff!(model::AbstractModel)

Compute the effective drag coefficient.
"""
function update_βeff!(model::AbstractModel)
    @unpack gh=model.fields
    gh.βeff .= gh.β ./ (1.0 .+ gh.quad_f2 .* gh.β)
    return model
end

"""
    update_βeff_on_uv_grids!(model::AbstractModel)

Interpolate the effective drag coefficient onto u- and v-grids, accounting for grounded fraction.
"""
function update_βeff_on_uv_grids!(model::AbstractModel)
    @unpack gh,gu,gv=model.fields
    @assert eltype(gh.grounded_fraction)==eltype(gh.βeff)

    T=eltype(gh.grounded_fraction)

    onesvec=ones(T,gh.nxh*gh.nyh)
    gu.βeff[gu.mask].=(gu.samp*(gu.cent'*(gh.crop*gh.βeff[:])))./(gu.samp*(gu.cent'*(gh.crop*onesvec)))
    ipolgfu=zeros(T,gu.nxu,gu.nyu);
    ipolgfu[gu.mask].=(gu.samp*(gu.cent'*(gh.crop*gh.grounded_fraction[:])))./(gu.samp*(gu.cent'*(gh.crop*onesvec)))
    gu.βeff[ipolgfu .> zero(T)] .= gu.βeff[ipolgfu .> zero(T)].*gu.grounded_fraction[ipolgfu .> zero(T)]./
                                                        ipolgfu[ipolgfu .> zero(T)]

    gv.βeff[gv.mask].=(gv.samp*(gv.cent'*(gh.crop*gh.βeff[:])))./(gv.samp*(gv.cent'*(gh.crop*onesvec)))
    ipolgfv=zeros(T,gv.nxv,gv.nyv);
    ipolgfv[gv.mask].=(gv.samp*(gv.cent'*(gh.crop*gh.grounded_fraction[:])))./(gv.samp*(gv.cent'*(gh.crop*onesvec)))
    gv.βeff[ipolgfv .> zero(T)] .= gv.βeff[ipolgfv .> zero(T)].*gv.grounded_fraction[ipolgfv .> zero(T)]./
                                                 ipolgfv[ipolgfv .> zero(T)];

    return model
end

"""
    update_rheological_operators!(model::AbstractModel)

Precompute various diagonal matrices used in defining the momentum operator.
"""
function update_rheological_operators!(model::AbstractModel)
    @unpack gh,gu,gv=model.fields
    @unpack params=model
    gh.dneghηav[] .= gh.crop*Diagonal(-gh.h[:].*gh.ηav[:])*gh.crop
    gu.dnegβeff[] .= gu.crop*Diagonal(-gu.βeff[:])*gu.crop
    gv.dnegβeff[] .= gv.crop*Diagonal(-gv.βeff[:])*gv.crop
    gh.dimplicit[] .= gh.crop*Diagonal(-params.density_ice * params.g * params.dt * gh.dsdh[:])*gh.crop
    return model
end

"""
    get_op(model::AbstractModel{T,N}) where {T,N}

 Get operator, defined as a LinearMap type.

"""
function get_op(model::AbstractModel{T,N}) where {T,N}
    @unpack gu,gv=model.fields
    n = gu.n + gv.n
    op_fun(x) = opvec(model,x)
    op=LinearMap{T}(op_fun,n;issymmetric=true,ismutating=false,ishermitian=true,isposdef=true)
end


text:
println("Hello wordl again!")

# something on main brach
# change on main branch

# somethig...

# and somenthing again to change it

max_stars_repo_name: Mousejr/testGit

max_issues_repo_name: Mousejr/testGit

max_forks_repo_name: Mousejr/testGit


text:
using CausalForest
using Test

@testset "oob regression random forest test" begin

    n, m = 10^3, 5
    features = randn(n, m)
    weights = rand(-2:2, m)
    labels = features * weights
    model = build_forest_oob(labels, features, -1, 100)
    predictions = apply_forest_oob(model)
    comparison = isapprox.(labels, predictions, atol=2.5)
    @test (sum(comparison)/n) > 0.8

end


text:
""" Shorthand for @artifact_str """
function artifact_helper(name::AbstractString)

    # The test here is whether to include the keyword quiet_download
    # The progress bar on artifacts >1.5 is wonky (especially when generating docs)
    if (VERSION > v"1.5")
        joinpath(ensure_artifact_installed(name, find_artifacts_toml(@__DIR__),
                                       quiet_download = true), name) * ".par"
    else
        joinpath(ensure_artifact_installed(name, find_artifacts_toml(@__DIR__)), name) * ".par"
    end

end

""" Given a molecule name and a database name, retrieve the transition file """
function artifact(molecule::AbstractString,
                  database::AbstractString = "hitran")
    @assert (database == "hitran") "This package currently only supports the hitran database!"
    @assert molecule in Absorption.mol_names "Not a supported molecule!"
    mol_id = Absorption.mol_number(molecule)
    return artifact_helper("hitran_molec_id_$(mol_id)_$(molecule)")
end


text:
@symbol_func function K_nT(cur_reactor::AbstractReactor)
  cur_K = 0.1581

  cur_K *= ( 1 + cur_reactor.f_D )

  cur_K *= ( 1 + cur_reactor.nu_n )

  cur_K *= ( 1 + cur_reactor.nu_T )

  cur_K /= ( 1 + cur_reactor.nu_n + cur_reactor.nu_T )

  cur_K
end


text:
cconvert(T::Type, x::VulkanStruct) = x
cconvert(T::Type{<:Ptr}, x::AbstractVector{<:VulkanStruct{false}}) = getproperty.(x, :vks)
cconvert(T::Type{<:Ptr}, x::AbstractVector{<:VulkanStruct{true}}) = (x, getproperty.(x, :vks))
cconvert(T::Type{<:Ptr}, x::VulkanStruct{false}) = Ref(x.vks)
cconvert(T::Type{<:Ptr}, x::VulkanStruct{true}) = (x, Ref(x.vks))
cconvert(T::Type{Ptr{_T}}, x::Vector{<:HighLevelStruct}) where {_T} = cconvert(T, convert(Vector{getproperty(@__MODULE__, Symbol(:_, nameof(eltype(x))))}, x))
cconvert(T::Type{Ptr{_T}}, x::HighLevelStruct) where {_T} = cconvert(T, convert(getproperty(@__MODULE__, Symbol(:_, nameof(typeof(x)))), x))
cconvert(T::Type{<:Ptr}, x::Handle) = x

convert(T::Type{<:Ptr}, x::Handle) = x.vks

unsafe_convert(T::Type, x::VulkanStruct) = x.vks
unsafe_convert(T::Type, x::Tuple{<:VulkanStruct{true}, <:Ref}) = unsafe_convert(T, last(x))
unsafe_convert(T::Type, x::Tuple{<:AbstractVector{<:VulkanStruct{true}}, <:Any}) = unsafe_convert(T, last(x))

"""
    `pointer_length(val)`

Return the length `val` considering it as an array.

Differ from `Base.length` in that `pointer_length(C_NULL) == 0` and that a `RefArray`s return the length of their roots.
"""
function pointer_length end

pointer_length(arr::Ptr{Nothing}) = 0
pointer_length(arr::AbstractArray) = length(arr)
pointer_length(arr::RefArray) = length(arr.roots)
pointer_length(arr::Tuple{<:Any,<:Any}) = length(first(arr))

convert_nonnull(T, val) = convert(T, val)
convert_nonnull(T, val::Ptr{Cvoid}) = val == C_NULL ? val : convert(T, val)

# # Initialization.

@nospecialize

initialize(::Type{T}) where {T<:VkCore.CEnum.Cenum} = typemin(T)
initialize(::Type{T}) where {T} = zero(T)
initialize(::Type{NTuple{N,T}}) where {N,T} = ntuple(Returns(zero(T)), N)
initialize(::Type{T}) where {T<:AbstractArray} = T()
initialize(::Type{<:Ptr}) = C_NULL

function init_empty(T)
    @assert isbitstype(T)
    ptr = Libc.calloc(1, sizeof(T))
    res = unsafe_load(Ptr{T}(ptr))
    Libc.free(ptr)
    res
end

function _initialize_core(T, refs)
    res = init_empty(T)
    sType = structure_type(T)
    if isempty(refs)
        pNext = C_NULL
    else
        ref = popfirst!(refs)
        ref[] = _initialize_core(Base.inferencebarrier(typeof(ref[])), refs)
        pNext = Base.unsafe_convert(Ptr{Cvoid}, ref)
    end
    setproperties(res, (; sType, pNext))
end

initialize_core(T, refs) = _initialize_core(Base.inferencebarrier(T), Base.inferencebarrier(refs))

# """
# Initialize a core structure of type `T` recursively, converting references provided in `vars_types` to pointers for their `next` chain.
# """
# function initialize_core(T::DataType, refs = [])
#     # This function takes a while to compile, if there are ways to speed it up that would be nice.
#     args = []
#     for (name, t) in zip(fieldnames(T), fieldtypes(T))
#         if name == :sType
#             push!(args, structure_type(T))
#         elseif name == :pNext
#             if isempty(refs)
#                 push!(args, C_NULL)
#             else
#                 ref, refs... = refs
#                 ref[] = initialize_core(eltype(ref), refs)
#                 push!(args, Base.unsafe_convert(Ptr{Cvoid}, ref))
#             end
#         elseif isstructtype(t) && !(t <: NTuple)
#             push!(args, initialize_core(t))
#         else
#             push!(args, initialize(t))
#         end
#     end
#     T(args...)
# end

"""
Initialize an intermediate structure, with a `next` chain built from initialized `Tnext` structs.
"""
function _initialize(T::Type{<:VulkanStruct}, Tnext)
    refs = Any[]
    for t in Tnext
        push!(refs, Ref{core_type(t)}())
    end
    vks = _initialize_core(Base.inferencebarrier(core_type(T)), Base.inferencebarrier(refs))
    T(vks, refs)
end

"""
Initialize a high-level structure, with a `next` chain built from initialized `Tnext` structs.
"""
function _initialize(T::Type{<:HighLevelStruct}, Tnext)
    args = []
    for (name, t) in zip(fieldnames(T), fieldtypes(T))
        if t === StructureType
            push!(args, structure_type(T))
        elseif name == :next
            push!(args, isempty(Tnext) ? C_NULL : _initialize(Base.inferencebarrier(Tnext[1]), Base.inferencebarrier(Tnext[2:end])))
        else
            push!(args, initialize(t))
        end
    end
    T(args...)
end

initialize(T::Type{<:HighLevelStruct}, args...) = _initialize(Base.inferencebarrier(T), Base.inferencebarrier(collect(args)))
initialize(T::Type{<:VulkanStruct}, args...) = _initialize(Base.inferencebarrier(T), Base.inferencebarrier(collect(args)))

# """
# Initialize a core structure of type `T` recursively, converting references provided in `vars_types` to pointers for their `next` chain.
# """
# function initialize_core(T::DataType, vars_types::Pair{Symbol,DataType}...)
#     # This function takes a while to compile, if there are ways to speed it up that would be nice.
#     args = []
#     for (name, t) in zip(fieldnames(T), fieldtypes(T))
#         if name == :sType
#             push!(args, structure_type(T))
#         elseif name == :pNext
#             if isempty(vars_types)
#                 push!(args, C_NULL)
#             else
#                 ((var, type), vars_types...) = vars_types
#                 ref_ex = initialize_core(type, (vars_types)...)
#                 push!(args, quote
#                     $var[] = $ref_ex
#                     Base.unsafe_convert(Ptr{Cvoid}, $var)
#                 end)
#             end
#         elseif isstructtype(t) && !(t <: NTuple)
#             push!(args, initialize_core(t))
#         else
#             push!(args, initialize(t))
#         end
#     end
#     :($T($(args...)))
# end

# """
# Initialize an intermediate structure, with a `next` chain built from initialized `Tnext` structs.
# """
# @generated function initialize(::Type{T}, Tnext...) where {T<:VulkanStruct}
#     Tnext = getindex.(getproperty.(Tnext, :parameters), 1)
#     vars_types = Symbol.(:ref_, 1:length(Tnext)) .=> core_type.(Tnext)
#     body = [:($var = Ref{$t}()) for (var, t) in vars_types]
#     vks = initialize_core(core_type(T), vars_types...)
#     Expr(:block, body..., :($T($vks, Any[$(first.(vars_types)...)])))
# end

# """
# Initialize a high-level structure, with a `next` chain built from initialized `Tnext` structs.
# """
# @generated function initialize(::Type{T}, Tnext...) where {T<:HighLevelStruct}
#     Tnext = getindex.(getproperty.(Tnext, :parameters), 1)
#     args = []
#     for (name, t) in zip(fieldnames(T), fieldtypes(T))
#         if t === StructureType
#             push!(args, structure_type(T))
#         elseif name == :next
#             push!(args, isempty(Tnext) ? C_NULL : initialize(Tnext...))
#         else
#             push!(args, initialize(t))
#         end
#     end
#     Expr(:call, nameof(T), args...)
# end

# @generated function load_next_chain(ptr, next_types...)
#     isempty(next_types) && return C_NULL
#     next_types = only.(getproperty.(next_types, :parameters))
#     T, Ts... = next_types
#     quote
#         @assert ptr ≠ C_NULL
#         vks = Base.unsafe_load(Base.unsafe_convert($(Ptr{core_type(T)}), ptr))
#         $T(vks, $(Ts...))
#     end
# end

@noinline function _load_next_chain(ptr, next_types)
    isempty(next_types) && return C_NULL
    T, Ts... = next_types
    @assert ptr ≠ C_NULL
    vks = Base.unsafe_load(Ptr{core_type(T)}(ptr))
    T(vks, Ts...)
end

load_next_chain(ptr, args...) = _load_next_chain(Base.inferencebarrier(ptr), Base.inferencebarrier(collect(args)))

@specialize

max_stars_repo_name: JuliaGPU/VulkanAbstraction

max_issues_repo_name: JuliaGPU/VulkanAbstraction

max_forks_repo_name: JuliaGPU/VulkanAbstraction


text:
module StructArrays

using Base: tuple_type_cons, tuple_type_head, tuple_type_tail, tail

import Requires
using PooledArrays: PooledArray

export StructArray, StructVector, LazyRow, LazyRows
export collect_structarray, fieldarrays
export replace_storage

include("interface.jl")
include("structarray.jl")
include("utils.jl")
include("collect.jl")
include("sort.jl")
include("groupjoin.jl")
include("lazy.jl")

function __init__()
    Requires.@require Tables="bd369af6-aec1-5ad0-b16a-f7cc5008161c" include("tables.jl")
    Requires.@require WeakRefStrings="ea10d353-3f73-51f8-a26c-33c1cb351aa5" begin
        function refs(a::WeakRefStrings.StringArray{T}) where {T}
            S = Union{WeakRefStrings.WeakRefString{UInt8}, typeintersect(T, Missing)}
            convert(WeakRefStrings.StringArray{S}, a)
        end
    end
end

end # module


text:
using Distributions
"""
    mutable struct MDP{T,R}
        observationspace::DiscreteSpace
        actionspace::DiscreteSpace
        state::Int64
        trans_probs::Array{T, 2}
        reward::R
        initialstates::Array{Int64, 1}
        isterminal::Array{Int64, 1}

A Markov Decision Process with `ns` states, `na` actions, current `state`,
`na`x`ns` - array of transition probabilites `trans_props` which consists for
every (action, state) pair of a (potentially sparse) array that sums to 1 (see
[`getprobvecrandom`](@ref), [`getprobvecuniform`](@ref),
[`getprobvecdeterministic`](@ref) for helpers to constract the transition
probabilities) `na`x`ns` - array of `reward`, array of initial states
`initialstates`, and `ns` - array of 0/1 indicating if a state is terminal.
"""
mutable struct MDP{T,R}
    observationspace::DiscreteSpace
    actionspace::DiscreteSpace
    state::Int64
    trans_probs::Array{T, 2}
    reward::R
    initialstates::Array{Int64, 1}
    isterminal::Array{Int64, 1}
end
function MDP(ospace, aspace, state, trans_probs::AbstractArray{T, 2},
             reward::R, initialstates, isterminal) where {T, R}
    if R <: AbstractMatrix
        reward = DeterministicStateActionReward(reward)
    end
    MDP{T,typeof(reward)}(ospace, aspace, state, trans_probs, reward, initialstates, isterminal)
end
function interact!(env::MDP, action)
    oldstate = env.state
    # @show oldstate
    run!(env, action)
    # @show env.state
    # @show env.isterminal[env.state]
    r = reward(env.reward, oldstate, action, env.state)
    (observation = env.state, reward = r, isdone = env.isterminal[env.state] == 1)
end
function getstate(env::MDP)
    (observation = env.state, isdone = env.isterminal[env.state] == 1)
end
function reset!(env::MDP)
    env.state = rand(ENV_RNG, env.initialstates)
    (observation = env.state, )
end

mutable struct ChangeMDP{TMDP}
    ns::Int64
    actionspace::DiscreteSpace
    changeprobability::Float64
    stochasticity::Float64
    mdp::TMDP
    switchflag::Array{Bool, 2}
    seed::Any
    rng::MersenneTwister # Used only for switches!
end
function ChangeMDP(; ns = 10, na = 4, nrewards = 2,
                    changeprobability = .01, stochasticity = 0.1,
                    seed = 3)
    rng = MersenneTwister(seed)
    mdpbase = MDP(ns, na, nrewards, init = "random")

    T = [rand(ENV_RNG, Dirichlet(ns, stochasticity)) for a in 1:na, s in 1:ns]
    T = removeautoconnections!(T, stochasticity, rng)
    T = replaceNaNswithdeterminism!(T, rng)
    T,~,~ = ensuretransitiontoterminal!(T, stochasticity, rng, mdpbase.isterminal)

    mdpbase.trans_probs = deepcopy(T)
    switchflag = Array{Bool, 2}(undef, na, ns)
    switchflag .= false
    ChangeMDP(ns, DiscreteSpace(na, 1), changeprobability, stochasticity,
                mdpbase, switchflag, seed, rng)
end
export ChangeMDP
getstate(env::ChangeMDP) = getstate(env.mdp)
reset!(env::ChangeMDP) = reset!(env.mdp)
function interact!(env::ChangeMDP, action)
    # # -----------------------------------------------------------
    # # # # For all s-a pairs, change them with prob pc: (option c.)
    # # -----------------------------------------------------------
    env.switchflag .= false
    r = rand!(env.rng, zeros(env.ns * env.mdp.actionspace.n))
    indicestoswitch = findall(r .< env.changeprobability)
    # @show r
    for i in indicestoswitch
        # @show i
        # println("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
        # println("%%%%%%%%%%%%%% CHANGE!!! %%%%%%%%%%%%%%%%%%%%%%")
        # println("%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%")
        # @show CartesianIndices(env.mdp.trans_probs)[i]
        # @show env.mdp.trans_probs[i]
        T = rand(env.rng, Dirichlet(env.ns, env.stochasticity))
        T = removeautoconnections!(T, env.stochasticity, env.rng, CartesianIndices(env.mdp.trans_probs)[i][2])
        T = replaceNaNswithdeterminism!(T, env.rng, CartesianIndices(env.mdp.trans_probs)[i][2])

        if ~all(env.mdp.trans_probs[i] .== T) # If there was indeed a change
            # @show CartesianIndices(env.mdp.trans_probs)[i]
            # @show env.mdp.trans_probs[i]
            env.mdp.trans_probs[i] = deepcopy(T)
            # @show env.mdp.trans_probs[i]
            env.switchflag[i] = true
        end
    end
    # Make sure there is at least one transtion to a terminalstate
    env.mdp.trans_probs, previousaction, previousstate = ensuretransitiontoterminal!(
                                                        env.mdp.trans_probs,
                                                        env.stochasticity,
                                                        env.rng,
                                                        env.mdp.isterminal)
    # @show previousaction, previousstate
    if previousaction!=0
        env.switchflag[previousaction, previousstate] = true
    end

    interact!(env.mdp, action)
end
# OLD versions:
# function interact!(env::ChangeMDP, action)
#     env.switchflag = false
#     r = rand(env.rng)
#     # @show r
#     if r < env.changeprobability # Switch or not!
#         # # -----------------------------------------------------------
#         # # # # Pick randomly an s-a pair and change it: (Problem: generative model is not the one we assume. Complicated)
#         # # -----------------------------------------------------------
#         a = rand(env.rng, 1:env.mdp.actionspace.n)
#         s = rand(env.rng, 1:env.ns)
#         T = rand(env.rng, Dirichlet(env.ns, env.stochasticity))
#         env.mdp.trans_probs[a, s] = deepcopy(T)
#         # # -----------------------------------------------------------
#         # # # # Change currect s-a pair: (This way env changes depend on agent's policy  -- We don't want this)
#         # # -----------------------------------------------------------
#         # env.mdp.trans_probs[action, env.mdp.state] = deepcopy(T)
#         # # -----------------------------------------------------------
#         # # # # Change the whole MDP: (Global changes -- We don't want this)
#         # # -----------------------------------------------------------
#         # T = [rand(env.rng, Dirichlet(env.ns, env.stochasticity))
#         #         for a in 1:env.mdp.actionspace.n, s in 1:env.ns]
#         # env.mdp.trans_probs = deepcopy(T)
#         env.switchflag = true
#     end
#     interact!(env.mdp, action)
# end
actionspace(env::ChangeMDP) = actionspace(env.mdp)
# ------------------------------------------------------------------------------
""" Jump MDP: similar to task used in fMRI project """
mutable struct JumpMDP{TMDP}
    ns::Int64
    actionspace::DiscreteSpace
    jumpprobability::Float64
    stochasticity::Float64
    mdp::TMDP
    switchflag::Array{Bool, 2} # It's rather a "jumpflag"
    seed::Any
    rng::MersenneTwister # Used only for switches!
    nonterminalstates::Array{Int, 1}
end
function JumpMDP(; ns = 10, na = 4, nrewards = 2,
                    jumpprobability = .01, stochasticity = 0.1,
                    seed = 3)
    rng = MersenneTwister(seed)
    mdpbase = MDP(ns, na, nrewards, init = "random")
    T = [rand(ENV_RNG, Dirichlet(ns, stochasticity)) for a in 1:na, s in 1:ns]
    mdpbase.trans_probs = deepcopy(T)
    switchflag = Array{Bool, 2}(undef, na, ns)
    switchflag .= false
    nonterminalstates = findall(mdpbase.isterminal .== 0)
    @show nonterminalstates
    JumpMDP(ns, DiscreteSpace(na, 1), jumpprobability, stochasticity,
                mdpbase, switchflag, seed, rng, nonterminalstates)
end
export JumpMDP
getstate(env::JumpMDP) = getstate(env.mdp)
reset!(env::JumpMDP) = reset!(env.mdp)
function interact!(env::JumpMDP, action)
    @show env.mdp.state
    @show action
    @show env.mdp.trans_probs[action, env.mdp.state]
    env.switchflag .= false
    r = rand(env.rng)
    if r < env.jumpprobability
        println("------ JUMP!")
        env.switchflag[action, env.mdp.state] = true
        interactjump!(env, action)
    else
        interact!(env.mdp, action)
    end
end
function interactjump!(env::JumpMDP, action)
    oldstate = env.mdp.state
    @show oldstate

    # possiblenextstates = deleteat!(collect(1:env.ns), collect(1:env.ns) .== oldstate)
    @show env.nonterminalstates
    possiblenextstates = copy(env.nonterminalstates)
    possiblenextstates = deleteat!(possiblenextstates, possiblenextstates .== oldstate)

    @show possiblenextstates
    env.mdp.state = rand(env.rng, possiblenextstates)
    @show env.mdp.state
    # @show env.isterminal[env.state]
    r = reward(env.mdp.reward, oldstate, action, env.mdp.state)
    @show r
    (observation = env.mdp.state, reward = r, isdone = env.mdp.isterminal[env.mdp.state] == 1)
end
actionspace(env::JumpMDP) = actionspace(env.mdp)
# ------------------------------------------------------------------------------
"""
    struct DeterministicNextStateReward
        value::Vector{Float64}
"""
struct DeterministicNextStateReward
    value::Vector{Float64}
end
reward(r::DeterministicNextStateReward, s, a, s′) = r.value[s′]
"""
    struct DeterministicStateActionReward
        value::Array{Float64, 2}

`value` should be a `na × ns`-matrix.
"""
struct DeterministicStateActionReward
    value::Array{Float64, 2}
end
reward(r::DeterministicStateActionReward, s, a, s′) = r.value[a, s]
"""
    struct NormalNextStateReward
        mean::Vector{Float64}
        std::Vector{Float64}
"""
struct NormalNextStateReward
    mean::Vector{Float64}
    std::Vector{Float64}
end
reward(r::NormalNextStateReward, s, a, s′) = r.mean[s′] + randn(ENV_RNG) * r.std[s′]
"""
    struct NormalStateActionReward
        mean::Array{Float64, 2}
        std::Array{Float64, 2}

`mean` and `std` should be `na × ns`-matrices.
"""
struct NormalStateActionReward
    mean::Array{Float64, 2}
    std::Array{Float64, 2}
end
reward(r::NormalStateActionReward, s, a, s′) = r.mean[a, s] + randn(ENV_RNG) * r.std[a, s]
"""
    getprobvecrandom(n)

Returns an array of length `n` that sums to 1. More precisely, the array is a
sample of a [Dirichlet
distribution](https://en.wikipedia.org/wiki/Dirichlet_distribution) with `n`
categories and ``α_1 = ⋯  = α_n = 1``.
"""
getprobvecrandom(n) = normalize(rand(ENV_RNG, n), 1)
"""
    getprobvecrandom(n, min, max)

Returns an array of length `n` that sums to 1 where all elements outside of
`min`:`max` are zero.
"""
getprobvecrandom(n, min, max) = sparsevec(collect(min:max),
                                          getprobvecrandom(max - min + 1), n)
"""
    getprobvecuniform(n)  = fill(1/n, n)
"""
getprobvecuniform(n) = fill(1/n, n)
"""
    getprobvecdeterministic(n, min = 1, max = n)

Returns a `SparseVector` of length `n` where one element in `min`:`max` has
value 1.
"""
getprobvecdeterministic(n, min = 1, max = n) = sparsevec([rand(ENV_RNG, min:max)], [1.], n)
# constructors
"""
    MDP(ns, na; init = "random")
    MDP(; ns = 10, na = 4, init = "random")

Return MDP with `init in ("random", "uniform", "deterministic")`, where the
keyword init determines how to construct the transition probabilites (see also
[`getprobvecrandom`](@ref), [`getprobvecuniform`](@ref),
[`getprobvecdeterministic`](@ref)).
"""
function MDP(ns, na; init = "random")
    r = randn(ENV_RNG, na, ns)
    func = eval(Meta.parse("getprobvec" * init))
    T = [func(ns) for a in 1:na, s in 1:ns]
    MDP(DiscreteSpace(ns, 1), DiscreteSpace(na, 1), rand(ENV_RNG, 1:ns), T, r,
        1:ns, zeros(ns))
end
MDP(; ns = 10, na = 4, init = "random") = MDP(ns, na, init = init)
function MDP(ns, na, r::Vector{Float64}; init = "random")
    func = eval(Meta.parse("getprobvec" * init))
    T = [func(ns) for a in 1:na, s in 1:ns]
    reward = DeterministicNextStateReward(r)
    MDP(DiscreteSpace(ns, 1), DiscreteSpace(na, 1), rand(ENV_RNG, 1:ns), T, reward,
        1:ns, zeros(ns))
end
function MDP(ns, na, nrewards::Int; init = "random")
    r = zeros(ns)
    [r[i] = 1. for i in rand(ENV_RNG, 1:ns, nrewards)]
    reward = DeterministicNextStateReward(r)
    isterminal = Int.(r)
    initialstates = 1:ns
    initialstates = @. initialstates[isterminal == 0] # Exclude terminalstates
    func = eval(Meta.parse("getprobvec" * init))
    T = [func(ns) for a in 1:na, s in 1:ns]
    MDP(DiscreteSpace(ns, 1), DiscreteSpace(na, 1), rand(ENV_RNG, 1:ns), T, reward,
        initialstates, isterminal)
end
actionspace(env::MDP) = env.actionspace

"""
    treeMDP(na, depth; init = "random", branchingfactor = 3)

Returns a tree structured MDP with na actions and `depth` of the tree.
If `init` is random, the `branchingfactor` determines how many possible states a
(action, state) pair has. If `init = "deterministic"` the `branchingfactor =
na`.
"""
function treeMDP(na, depth;
                 init = "random",
                 branchingfactor = 3)
    isdet = (init == "deterministic")
    if isdet
        branchingfactor = na
        ns = na.^(0:depth - 1)
    else
        ns = branchingfactor.^(0:depth - 1)
    end
    cns = cumsum(ns)
    func = eval(Meta.parse("getprobvec" * init))
    T = Array{SparseVector, 2}(undef, na, cns[end])
    for i in 1:depth - 1
        for s in 1:ns[i]
            for a in 1:na
                lb = cns[i] + (s - 1) * branchingfactor + (a - 1) * isdet + 1
                ub = isdet ? lb : lb + branchingfactor - 1
                T[a, (i == 1 ? 0 : cns[i-1]) + s] = func(cns[end] + 1, lb, ub)
            end
        end
    end
    r = zeros(na, cns[end] + 1)
    isterminal = [zeros(cns[end]); 1]
    for s in cns[end-1]+1:cns[end]
        for a in 1:na
            r[a, s] = -rand(ENV_RNG)
            T[a, s] = getprobvecdeterministic(cns[end] + 1, cns[end] + 1,
                                              cns[end] + 1)
        end
    end
    MDP(DiscreteSpace(cns[end] + 1, 1), DiscreteSpace(na, 1), 1, T, r, 1:1, isterminal)
end

function emptytransprob!(v::SparseVector)
    empty!(v.nzind); empty!(v.nzval)
end
emptytransprob!(v::Array{Float64, 1}) = v[:] .*= 0.

"""
    setterminalstates!(mdp, range)

Sets `mdp.isterminal[range] .= 1`, empties the table of transition probabilities
for terminal states and sets the reward for all actions in the terminal state to
the same value.
"""
function setterminalstates!(mdp, range)
    mdp.isterminal[range] .= 1
    for s in findall(mdp.isterminal)
        mdp.reward[:, s] .= mean(mdp.reward[:, s])
        for a in 1:mdp.na
            emptytransprob!(mdp.trans_probs[a, s])
        end
    end
end

# run MDP

"""
    run!(mdp::MDP, action::Int64)

Transition to a new state given `action`. Returns the new state.
"""
function run!(mdp::MDP, action::Int64)
    if mdp.isterminal[mdp.state] == 1
        reset!(mdp)
    else
        mdp.state = wsample(ENV_RNG, mdp.trans_probs[action, mdp.state])
        (observation = mdp.state,)
    end
end

"""
    run!(mdp::MDP, policy::Array{Int64, 1}) = run!(mdp, policy[mdp.state])

"""
run!(mdp::MDP, policy::Array{Int64, 1}) = run!(mdp, policy[mdp.state])

function replaceNaNswithdeterminism!(T, rng)
    nanvectorindices = findall([any(isnan.(i)) for i in T])
    if ~isempty(nanvectorindices)
        for k in nanvectorindices
            # @show k
            T[k[1], k[2]]  = replaceNaNswithdeterminism!(T[k[1], k[2]], rng, k[2])
        end
    end
    T
end
function replaceNaNswithdeterminism!(T::Array{Float64,1}, rng, stateindex)
    ns = length(T)
    # @show T
    if any(isnan.(T))
        # @show stateindex
        T = zeros(ns)
        possiblenextstates = deleteat!(collect(1:ns), stateindex)
        nextstate = rand(rng, possiblenextstates)
        T[nextstate] = 1.
    end
    T
end
export replaceNaNswithdeterminism!
function removeautoconnections!(T, stochasticity, rng)
    # na = size(T,1); ns = size(T,2)
    for a in 1:size(T,1)
        for s in 1:size(T,2)
            # @show (a, s)
            T[a, s] = removeautoconnections!(T[a, s], stochasticity, rng, s)
        end
    end
    T
end
function removeautoconnections!(T::Array{Float64,1}, stochasticity, rng, stateindex)
    ns = length(T)
    # @show T
    while T[stateindex] > 0.9999
        # @show stateindex
        T = rand(rng, Dirichlet(ns, stochasticity))
    end
    T
end
export removeautoconnections!
function ensuretransitiontoterminal!(T, stochasticity, rng, isterminal)
    previousaction = 0; previousstate = 0
    terminalstates = findall(y-> y==1, isterminal)
    istheretransition = checktransitiontoterminal(T, terminalstates)
    # @show typeof(istheretransition)
    if ~istheretransition
        # println("i m here")
        possibleprevioustates = deleteat!(collect(1:size(T,2)), terminalstates)
        previousstate = rand(rng, possibleprevioustates)
        previousaction = rand(rng, collect(1:size(T,1)))
        chosenterminalstate = rand(rng, terminalstates)
        Ttemp = zeros(size(T,2))
        while Ttemp[chosenterminalstate] < 0.1
            # println("i m here too")
            # @show stateindex
            Ttemp = rand(rng, Dirichlet(size(T,2), stochasticity))
            if any(isnan.(Ttemp))
                Ttemp = zeros(size(T,2))
                Ttemp[chosenterminalstate] = 1.
            end
        end
        T[previousaction, previousstate] = copy(Ttemp)
    end
    T, previousaction, previousstate
end
export ensuretransitiontoterminal!
function checktransitiontoterminal(T, terminalstates)
    istheretransition = false
    # @show typeof(istheretransition)
    for i in terminalstates
        for a in 1:size(T,1), s in 1:size(T,2)
            # @show a
            # @show s
            # @show T[a, s][i]
            if T[a, s][i] > 0.1
                istheretransition = true
                break;
            end
            # @show istheretransition
            # @show typeof(istheretransition)
        end
        if istheretransition
            break;
        end
    end
    istheretransition
end
export checktransitiontoterminal


text:
## CONSTANTS

const Dist = Distributions

## EQUALITY OF DISTRIBUTIONS (minor type piracy)

function ==(d1::D, d2::D) where D<:Dist.Sampleable
    ret = true
    for fld in fieldnames(D)
        ret = ret && getfield(d1, fld) == getfield(d2, fld)
    end
    return ret
end

# DISTRIBUTION AS TRAIT (needed?)

# fallback:
isdistribution(d) = isdistribution(typeof(d))
isdistribution(::Type{<:Any}) = false

# for anything sampleable in Distributions.jl:
isdistribution(::Type{<:Dist.Sampleable}) = true

# ADD TO distributions.jl TYPE HIERARCHY TO ACCOUNT FOR NON-EUCLIDEAN
# SUPPORTS

abstract type NonEuclidean <: Distributions.ValueSupport end

## UNIVARIATE NOMINAL PROBABILITY DISTRIBUTION

# to pad a dictionary of probabilities with zeros for unseen levels,
# if the key type is CategoricalElement:
pad_probabilities(prob_given_level) = prob_given_level
function pad_probabilities(prob_given_level::Dict{<:CategoricalElement, T}) where T
    proto_element = first(keys(prob_given_level))
    return merge!(Dict([y=>zero(T) for y in classes(proto_element)]),
                  prob_given_level)
end

"""
    UnivariateFinite(levels, p)

A discrete univariate distribution whose finite support is the
elements of the vector `levels`, and whose corresponding probabilities
are elements of the vector `p`, which must sum to one.

In the special case that `levels` has type `AbstractVector{L}` where
`L <: CategoricalValue` or `L <: CategoricalString` (for example
`levels` is a `CategoricalVector`) the constructor adds the unobserved
classes (from the common pool) with probability zero.

    UnivariateFinite(prob_given_level)

A discrete univariate distribution whose finite support is the set of
keys of the provided dictionary, `prob_given_level`. The dictionary
values specify the corresponding probabilities, which must be
nonnegative and sum to one.

In the special case that `keys(prob_given_level)` has type
`AbstractVector{L}` where `L <: CategoricalValue` or `L <:
CategoricalString` (for example it is a `CategoricalVector`) the
constructor adds the unobserved classes from the common pool with
probability zero.

    levels(d::UnivariateFinite)

Return the levels of `d`.

````julia
d = UnivariateFinite(["yes", "no", "maybe"], [0.1, 0.2, 0.7])
pdf(d, "no") # 0.2
mode(d) # "maybe"
rand(d, 5) # ["maybe", "no", "maybe", "maybe", "no"]
d = fit(UnivariateFinite, ["maybe", "no", "maybe", "yes"])
pdf(d, "maybe") ≈ 0.5 # true
levels(d) # ["yes", "no", "maybe"]
````

If the element type of `v` is a `CategoricalValue` or
`CategoricalString`, then `fit(UnivariateFinite, v)` assigns a
probability of zero to unobserved classes from the common pool.

See also `classes`.

"""
struct UnivariateFinite{L,T<:Real} <: Dist.Distribution{Dist.Univariate,NonEuclidean}
    prob_given_level::Dict{L,T}
    function UnivariateFinite{L,T}(prob_given_level::Dict{L,T}) where {L,T<:Real}
        p = values(prob_given_level) |> collect
        Dist.@check_args(UnivariateFinite, Dist.isprobvec(p))
        return new{L,T}(pad_probabilities(prob_given_level))
    end
end

UnivariateFinite(prob_given_level::Dict{L,T}) where {L,T<:Real} =
    UnivariateFinite{L,T}(prob_given_level)

function UnivariateFinite(levels::AbstractVector, p::AbstractVector{<:Real})
        Dist.@check_args(UnivariateFinite, length(levels)==length(p))
    prob_given_level = Dict([levels[i]=>p[i] for i in eachindex(p)])
    return  UnivariateFinite(prob_given_level)
end

CategoricalArrays.levels(d::UnivariateFinite) = collect(keys(d.prob_given_level))

function average(dvec::Vector{UnivariateFinite{L,T}}; weights=nothing) where {L,T}

    n = length(dvec)

    Dist.@check_args(UnivariateFinite, weights == nothing || n==length(weights))

    if weights == nothing
        weights = fill(1/n, n)
    else
        weights = weights/sum(weights)
    end

    # get all levels:
    levels = reduce(union, [keys(d.prob_given_level) for d in dvec])

    z = Dict{L,T}([x => zero(T) for x in levels]...)
    prob_given_level_vec = map(dvec) do d
        merge(z, d.prob_given_level)
    end

    # initialize the prob dictionary for the distribution sum:
    prob_given_level = Dict{L,T}()
    for x in levels
        prob_given_level[x] = zero(T)
    end

    # sum up:
    for x in levels
        for k in 1:n
            prob_given_level[x] += weights[k]*prob_given_level_vec[k][x]
        end
    end

    return UnivariateFinite(prob_given_level)

end

function Distributions.mode(d::UnivariateFinite)
    dic = d.prob_given_level
    p = values(dic)
    max_prob = maximum(p)
    m = first(first(dic)) # mode, just some level for now
    for (x, prob) in dic
        if prob == max_prob
            m = x
            break
        end
    end
    return m
end

Distributions.pdf(d::UnivariateFinite, x) = d.prob_given_level[x]

"""
    _cummulative(d::UnivariateFinite)

Return the cummulative probability vector `[0, ..., 1]` for the
distribution `d`, using whatever ordering is used in the dictionary
`d.prob_given_level`. Used only for to implement random sampling from
`d`.

"""
function _cummulative(d::UnivariateFinite{L,T}) where {L,T<:Real}
    p = collect(values(d.prob_given_level))
    K = length(p)
    p_cummulative = Array{T}(undef, K + 1)
    p_cummulative[1] = zero(T)
    p_cummulative[K + 1] = one(T)
    for i in 2:K
        p_cummulative[i] = p_cummulative[i-1] + p[i-1]
    end
    return p_cummulative
end

"""
_rand(p_cummulative)

Randomly sample the distribution with discrete support `1:n` which has
cummulative probability vector `p_cummulative=[0, ..., 1]` (of length
`n+1`). Does not check the first and last elements of `p_cummulative`
but does not use them either.

"""
function _rand(p_cummulative)
    real_sample = rand()
    K = length(p_cummulative)
    index = K
    for i in 2:K
        if real_sample < p_cummulative[i]
            index = i - 1
            break
        end
    end
    return index
end

function Base.rand(d::UnivariateFinite)
    p_cummulative = _cummulative(d)
    levels = collect(keys(d.prob_given_level))
    return levels[_rand(p_cummulative)]
end

function Base.rand(d::UnivariateFinite, n::Int)
    p_cummulative = _cummulative(d)
    levels = collect(keys(d.prob_given_level))
    return [levels[_rand(p_cummulative)] for i in 1:n]
end

function Distributions.fit(d::Type{<:UnivariateFinite}, v::AbstractVector)
    vpure = skipmissing(v) |> collect
    isempty(vpure) && error("No non-missing data to fit. ")
    N = length(vpure)
    count_given_level = Dist.countmap(vpure)
    prob_given_level = Dict([x=>c/N for (x, c) in count_given_level])
    return UnivariateFinite(prob_given_level)
end

## NORMAL DISTRIBUTION ARITHMETIC

# *(lambda::Number, d::Dist.Normal) = Dist.Normal(lambda*d.μ, lambda*d.σ)
# function +(ds::Dist.Normal...)
#     μ = sum([d.μ for d in ds])
#     σ = sum([d.σ^2 for d in ds]) |> sqrt
#     return Dist.Normal(μ, σ)
# end

# hack for issue #809 in Dist.jl:
Distributions.fit(::Type{Dist.Normal{T}}, args...) where T =
    Dist.fit(Dist.Normal, args...)


text:
module ParallelCoverallsSetup

export i_am_partially_tested, i_am_fully_tested, i_am_not_tested, i_will_throw

function i_am_partially_tested()
  local val::Int

  if Sys.islinux()
    val = 1
  elseif Sys.isapple()
    val = 2
  elseif Sys.iswindows()
    val = 3
  else
    val = 4
  end

  return val
end

function i_am_fully_tested()
  val = rand()
  return val
end

function i_am_not_tested(x)
  if x == 1
    return "x is one"
  end

  return "x is not one"
end

function i_will_throw()
  x = rand(10)
  mat = rand(10, 10)
  y = mat * x
  @assert y === x
end

end # module


text:
# Alias of Model{K,Kₑ,Kₘ, Kₑₜ,Kₘₜ, Kₑₗ,Kₘₗ, K², K₊₂,AK₊₂} for the waveguide mode
# equations for full 2D cross section.
#
# The type parameters AK₊₁ and AK₊₂ specify device-specific arrays types (e.g., CuArray) and
# are user-defined in the constructor.
const ModelFull{AK₊₂} = Model{2,3,3, 2,2, 1,1, 4, 4,AK₊₂}

# Convenience constructor
function ModelFull(grid::Grid; Atype::Type=Array)
    ft_eq = EE

    cmpₛ = SInt(1,2)  # shapes in transverse dimension (xy-plane)

    cmpₑₜ = SInt(1,2)  # transverse E-field (x- and y-components)
    cmpₘₜ = SInt(1,2)  # transverse H-field (x- and y-components)

    cmpₑₗ = SInt(3)  # longitudinal E-field (z-component)
    cmpₘₗ = SInt(3)  # longitudinal H-field (z-component)

    iseₜ˔shp = false
    ishₜ˔shp = false

    return ModelFull{Atype{ComplexF,4}}(;ft_eq, grid, cmpₛ, cmpₑₜ, cmpₘₜ, cmpₑₗ, cmpₘₗ, iseₜ˔shp, ishₜ˔shp)
end

# Assign the material parameters on the grid and smooth them.
function calc_matparams!(mdl::ModelFull)
    # Take necessary fields of mdl.
    grid = mdl.grid

    isbloch = grid.isbloch
    N = grid.N
    l = grid.l
    σ = grid.σ

    lg = grid.ghosted.l
    τl = grid.ghosted.τl
    ∆τ = grid.ghosted.∆τ

    boundft = mdl.boundft
    gtₑ = ft2gt.(EE,boundft)
    gtₘ = ft2gt.(HH,boundft)

    εₜarr = mdl.εₜarr
    µₜarr = mdl.μₜarr
    εₗarr = mdl.εₗarr
    µₗarr = mdl.μₗarr

    oind2shp = mdl.oind2shp
    oind2εind = mdl.oind2εind
    oind2μind = mdl.oind2μind
    εind2ε = mdl.εind2ε
    μind2μ = mdl.μind2μ

    cmpₑₜ = mdl.cmpₑₜ
    cmpₘₜ = mdl.cmpₘₜ
    cmpₑₗ = mdl.cmpₑₗ
    cmpₘₗ = mdl.cmpₘₗ

    εind2εₜ = sub_pind2matprm(εind2ε, cmpₑₜ)
    εind2εₗ = sub_pind2matprm(εind2ε, cmpₑₗ)
    μind2μₜ = sub_pind2matprm(μind2μ, cmpₘₜ)
    μind2μₗ = sub_pind2matprm(μind2μ, cmpₘₗ)

    iseₜ˔shp = mdl.iseₜ˔shp
    ishₜ˔shp = mdl.ishₜ˔shp
    iseₗ˔shp = true
    ishₗ˔shp = true

    # Create temporary storages.
    εxx_oind2d = create_oind_array(N)
    εyy_oind2d = create_oind_array(N)
    εzz_oind2d = create_oind_array(N)
    εoo_oind2d = create_oind_array(N)

    μxx_oind2d = create_oind_array(N)
    μyy_oind2d = create_oind_array(N)
    μzz_oind2d = create_oind_array(N)
    μoo_oind2d = create_oind_array(N)

    # Assign the material parameters.
    # (See the figure in L10 - Eigenmode Analysis > 1.3 Waveguide mode analysis > Matrix
    # equation formulation.)
    assign_param!(εₜarr, (εyy_oind2d,εxx_oind2d), oind2shp, oind2εind, εind2εₜ, gtₑ, τl, isbloch)  # diagonal entries of εₜ tensors
    assign_param!(εₜarr, tuple(μzz_oind2d), oind2shp, oind2εind, εind2εₜ, gtₑ, τl, isbloch)  # off-diagonal entries of εₜ tensors
    assign_param!(µₗarr, tuple(εoo_oind2d), oind2shp, oind2μind, μind2μₗ, gtₘ, τl, isbloch)  # μₗ tensors (rank-0, so scalars)

    assign_param!(µₜarr, (μyy_oind2d,μxx_oind2d), oind2shp, oind2μind, μind2μₜ, gtₘ, τl, isbloch)  # diagonal entries of μₜ tensors
    assign_param!(µₜarr, tuple(εzz_oind2d), oind2shp, oind2μind, μind2μₜ, gtₘ, τl, isbloch)  # off-diagonal entries of μₜ tensors
    assign_param!(εₗarr, tuple(μoo_oind2d), oind2shp, oind2εind, εind2εₗ, gtₑ, τl, isbloch)  # εₗ tensors (rank-0, so scalars)

    # Some temporary arrays should be identical.  For example, εyy_oind2d and μxx_oind2d are
    # evaluated at the locations surrounding the εyy and μxx locations, which are the same
    # locations.  (See the figure in L10 - Eigenmode Analysis > 1.3 Waveguide mode analysis
    # > Matrix equation formulation.)  Therefore, εyy_oind2d and μxx_oind2d are evaluated at
    # the same locations, so they should be the same object indices.
    @assert isequal(εxx_oind2d, μyy_oind2d)
    @assert isequal(εyy_oind2d, μxx_oind2d)
    @assert isequal(εoo_oind2d, εzz_oind2d)
    @assert isequal(μoo_oind2d, μzz_oind2d)

    # Smooth the material parameters.
    smooth_param!(εₜarr, (εxx_oind2d,εyy_oind2d), oind2shp, oind2εind, εind2εₜ, gtₑ, l, lg, σ, ∆τ, iseₜ˔shp)  # diagonal entries of εₜ tensors
    smooth_param!(εₜarr, tuple(εoo_oind2d), oind2shp, oind2εind, εind2εₜ, gtₑ, l, lg, σ, ∆τ, iseₜ˔shp)  # off-diagonal entries of εₜ tensors
    smooth_param!(µₗarr, tuple(μzz_oind2d), oind2shp, oind2μind, μind2μₗ, gtₘ, l, lg, σ, ∆τ, ishₗ˔shp)  # μₗ tensors (rank-0, so scalars)

    smooth_param!(µₜarr, (μxx_oind2d,μyy_oind2d), oind2shp, oind2μind, μind2μₜ, gtₘ, l, lg, σ, ∆τ, ishₜ˔shp)  # diagonal entries of μₜ tensors
    smooth_param!(µₜarr, tuple(μoo_oind2d), oind2shp, oind2μind, μind2μₜ, gtₘ, l, lg, σ, ∆τ, ishₜ˔shp)  # off-diagonal entries of μₜ tensors
    smooth_param!(εₗarr, tuple(εzz_oind2d), oind2shp, oind2εind, εind2εₗ, gtₑ, l, lg, σ, ∆τ, iseₗ˔shp)  # εₗ tensors (rank-0, so scalars)

    return nothing
end

function complete_fields(β::Number,  # propagation constant
                         fₜ::AbsVecNumber,  # calculated transverse field Fₜ as column vector
                         ft::FieldType,  # type of input field Fₜ
                         ω::Number,
                         Ps::Tuple22{AbsMatNumber},
                         ∇̽s::Tuple22{AbsMatNumber},
                         πcmps::Tuple2{AbsMatNumber},
                         mdl::ModelFull)
    ft′ = alter(ft)

    # Calculate the transverse components of the complementary field.
    βF′ₜgen = create_βF′ₜgen(ft, ω, Ps, ∇̽s, πcmps)
    f′ₜ = (βF′ₜgen * fₜ) ./ β

    # Calculate the longitudinal components of the main and complementary fields.
    iF′ₗgen = create_iF′ₗgen(ft, ω, Ps, ∇̽s)
    iFₗgen = create_iF′ₗgen(ft′, ω, Ps, ∇̽s)

    fₗ = (iFₗgen * f′ₜ) ./ im
    f′ₗ = (iF′ₗgen * fₜ) ./ im

    # Take the Cartesian components of the fields and reshape them into arrays.
    sz = size(mdl, ft)
    sz′ = size(mdl, ft′)
    @assert sz′ == sz

    fₜ = reshape(fₜ, sz)
    f′ₜ = reshape(f′ₜ, sz)

    order_cmpfirst = mdl.order_cmpfirst
    d = order_cmpfirst ? 1 : length(sz)

    Fx = selectdim(fₜ, d, 1)
    Fy = selectdim(fₜ, d, 2)

    F′x = selectdim(f′ₜ, d, 1)
    F′y = selectdim(f′ₜ, d, 2)

    sz_grid = size(mdl.grid)
    Fz = reshape(fₗ, sz_grid)
    F′z = reshape(f′ₗ, sz_grid)

    F = (Fx, Fy, Fz)
    F′ = (F′x, F′y, F′z)

    E, H = ft==EE ? (F,F′) : (F′,F)

    return (E=E, H=H)
end

# See the field configurations in L10 - Eigenmode Analysis > Matrix equation formulation.
function poynting(E::Tuple3{AbsMatNumber},  # (Ex, Ey, Ez)
                  H::Tuple3{AbsMatNumber},  # (Hx, Hy, Hz)
                  mdl::ModelFull)
    boundft = mdl.boundft
    grid = mdl.grid
    isbloch = grid.isbloch

    ∆l, ∆l⁻¹ = create_∆ls(grid, boundft)

    Ex, Ey, Ez = E
    Hx, Hy, Hz = H

    nX, nY = 1, 2

    Ftemp = ArrComplexF(undef, size(grid))  # matrix; temporary storage

    # Calculate Sx.
    Sx = ArrComplexF(undef, size(grid))  # matrix

    m̂xHz = Ftemp; interp_field!(m̂xHz, Hz, HH, nX, ∆l, ∆l⁻¹, boundft, isbloch)  # intepolate Hz at Ey-locations
    Sx .= Ey .* conj.(m̂xHz)  # EyHz evaluated at centers of x-normal edges

    m̂xHy = Ftemp; interp_field!(m̂xHy, Hy, HH, nX, ∆l, ∆l⁻¹, boundft, isbloch)  # intepolate Hy at Ez-locations
    m̂xHy .= Ez .* conj.(m̂xHy)  # EzHy evaluated at end points of x-normal edges; need to interpolate at centers
    Sx₋ = m̂xHy

    interp_field!(Sx, Sx₋, EE, nY, ∆l, ∆l⁻¹, boundft, isbloch, α=-1.0, overwrite=false)  # interpolate Sx₋ at centers of x-normal edges

    # Calculate Sy.
    Sy = ArrComplexF(undef, size(grid))  # matrix

    m̂yHz = Ftemp; interp_field!(m̂yHz, Hz, HH, nY, ∆l, ∆l⁻¹, boundft, isbloch)  # intepolate Hz at Ex-locations
    Sy .= -1.0 .* Ex .* conj.(m̂yHz)  # ExHz evaluated at centers of y-normal edges

    m̂yHx = Ftemp; interp_field!(m̂yHx, Hx, HH, nY, ∆l, ∆l⁻¹, boundft, isbloch)  # intepolate Hx at Ez-locations
    m̂yHx .= Ez .* conj.(m̂yHx)  # EzHx evaluated at end points of y-normal edges; need to interpolate at centers
    Sy₊ = m̂yHx

    interp_field!(Sy, Sy₊, EE, nX, ∆l, ∆l⁻¹, boundft, isbloch, overwrite=false)  # interpolate Sy₊ at centers of y-normal edges

    # Calculate Sz.
    Sz = ArrComplexF(undef, size(grid))  # matrix

    Sz₊ = Ftemp; Sz₊ .= Ex .* conj.(Hy)  # Ex and Hy are already co-located
    interp_field!(Sz, Sz₊, EE, nY, ∆l, ∆l⁻¹, boundft, isbloch)  # Sz₊ interpolated at centers of pixels

    Sz₋ = Ftemp; Sz₋ .= Ey .* conj.(Hx)  # Ey and Hx are already co-located
    interp_field!(Sz, Sz₋, EE, nX, ∆l, ∆l⁻¹, boundft, isbloch, α=-1.0, overwrite=false)  # Sz₋ interpolated at centers of pixels

    return Sx, Sy, Sz
end


text:
module ShaleDrillingLikelihood_Parameters

using Revise

DOBTIME = false
DOPROFILE = false

using ShaleDrillingLikelihood
using ShaleDrillingLikelihood.SDLParameters

using Test
using SparseArrays
using Random
using Distributions
using Dates
using StatsBase
using LinearAlgebra

using Base.Iterators: product, OneTo

using ShaleDrillingLikelihood: ObservationDrill,
    SimulationDraw,
    zero_out_small_probs,
    lrvar,
    lrstd,
    lrmean,
    simulate,
    theta_royalty_ρ,
    revenue,
    cost,
    extend,
    theta_ρ,
    vw_revenue,
    vw_cost,
    vw_extend,
    DataDynamicDrill,
    max_states,
    simloglik!,
    SimulationDraws,
    _model,
    theta_drill_ρ,
    update!,
    value_function,
    EV, dEV,
    _x, _y,
    idx_produce_ψ, idx_drill_ψ,
    idx_produce_g, idx_drill_g,
    merge_thetas,
    split_thetas,
    theta_produce,
    serial_simloglik!,
    parallel_simloglik!,
    drill,
    total_wells_drilled,
    _D

using ShaleDrillingLikelihood.SDLParameters: TestPriceDates,
    TestPriceGrid, TestPriceTransition, TestAR1Process,
    grid, series, transition

println("print to keep from blowing up")

@testset "Parameters" begin
    rwrd_c = DrillReward(DrillingRevenue(Constrained(;log_ogip=AlphaG(), α_ψ=AlphaPsi()),NoTrend(),NoTaxes()), DrillingCost_constant(), ExtensionCost_Constant())
    rwrd_u = UnconstrainedProblem(rwrd_c)
    rwrd = rwrd_u

    θ_royalty = Theta(RoyaltyModel())
    θ_produce = Theta(ProductionModel())
    θ_drill_c = Theta(rwrd_c)
    θ_drill_u = Theta(rwrd_u)

    @test _nparm(rwrd) == length(θ_drill_u)

    num_zt = 150          # drilling
    nz = 15
    num_i = 50
    nψ =  13
    M = 50

    # prices
    tpd = TestPriceDates(num_zt)
    proc = TestAR1Process()
    simz = simulate(proc, num_zt)
    etv = ExogTimeVars(tuple.(simz), tpd)
    grd = TestPriceGrid(simz, nz; process=proc)
    trans = TestPriceTransition(grd...; process=proc)
    PriceProcess(proc, grd, etv, trans)

    @test minimum(etv) isa Tuple
    @test length(product(grd...)) == nz

    price_process = TestPriceProcess(;process=proc, nsim=num_zt, ngrid=nz)
    _zchars = SDLParameters.zchars(price_process)

    # ogips
    ogip_dist = NormalOGIP()
    @test ogip_dist == Normal(4.68, 0.31)
    logOGIP(num_i, Normal())
    log_ogip = logOGIP(num_i)

    # observations
    obs_per_well = 10:20  # pdxn

    # roylaty
    royalty_rates = [1/8, 3/16, 1/4,]
    @test royalty_rates == RoyaltyRates(3)
    @test royalty_rates == RoyaltyRates(:small)
    @test RoyaltyRates(6) == RoyaltyRates(:full)
    @test length(RoyaltyRates(6)) == 6
    @test issorted(RoyaltyRates(6))
    @test issorted(RoyaltyRates(3))

    Theta(RoyaltyModel())

    RoyaltyExogVars(num_i, royalty_rates, θ_royalty)
    RoyaltyExogVars(num_i, royalty_rates, θ_royalty, log_ogip)
    Xroyalty = RoyaltyExogVars(num_i, royalty_rates, θ_royalty, log_ogip, RoyaltyModel())

    # geology
    ψs = PsiSpace(nψ)
    @test ψs == range(-4.5; stop=4.5, length=nψ)

    wp = FullLeasedProblem()

    u,v = randn(num_i), randn(num_i)

    # construct royalty data
    data_roy = DataRoyalty(u, v, Xroyalty, θ_royalty, royalty_rates)
    royrates_endog = observed_choices(data_roy)
    @test royalty_rates[_y(data_roy)] == royrates_endog
    royrates_exog = sample(royalty_rates, num_i)

    _ichars = SDLParameters.ichars(log_ogip, royrates_endog)
    @test _ichars == [gr for gr in zip(log_ogip, royrates_endog)]

    # model
    ddm_no_t1ev     = TestDynamicDrillModel(price_process, false; reward=rwrd_u)
    ddm_with_t1ev   = TestDynamicDrillModel(price_process, true;  reward=rwrd_u)
    ddm_c_no_t1ev   = TestDynamicDrillModel(price_process, false; reward=rwrd_c)
    ddm_c_with_t1ev = TestDynamicDrillModel(price_process, true;  reward=rwrd_c)

    # construct drilling data
    ddm_opts = (minmaxleases=1:1, nper_initial=40:40, tstart=1:50)
    data_drill_w = DataDynamicDrill(u, v, _zchars, _ichars, ddm_with_t1ev, θ_drill_u; ddm_opts...)
    data_drill_n = DataDynamicDrill(u, v, _zchars, _ichars, ddm_no_t1ev,   θ_drill_u; ddm_opts...)

    # constrained versions of above
    data_drill_w_con = DataDrill(data_drill_w, ddm_c_with_t1ev)
    data_drill_n_con = DataDrill(data_drill_n, ddm_c_no_t1ev)

    # number of wells drilled
    nwells_w = total_wells_drilled(data_drill_w)
    nwells_n = total_wells_drilled(data_drill_n)
    @test nwells_w == map(s -> _D(wp,s), max_states(data_drill_w))
    @test nwells_n == map(s -> _D(wp,s), max_states(data_drill_n))

    data_produce_w = DataProduce(u, nwells_w, obs_per_well, θ_produce, log_ogip)
    data_produce_n = DataProduce(u, nwells_n, obs_per_well, θ_produce, log_ogip)

    data_d_sm = data_drill_w_con
    data_d_lg = data_drill_w
    data_p    = data_produce_w
    data_r    = data_roy
    empty = EmptyDataSet()

    coef_links = CoefLinks(ddm_no_t1ev)
    @test coef_links == [(idx_produce_ψ, idx_drill_ψ,), (idx_produce_g, idx_drill_g)]

    data_dril = DataSetofSets(data_d_sm, empty, empty)
    data_full = DataSetofSets(data_d_lg, data_r, data_p, CoefLinks(data_d_lg))

    theta_dril = θ_drill_c
    theta_full = merge_thetas((θ_drill_u, θ_royalty, θ_produce), data_full)

end # test
end # module


text:
include("/home/jefter66/Workspace/lablib/src/calcs/Uncertainty.jl")

# Dados para calculo pela conservação de energia

# converti para o SI

# a massa é a soma das massas de cada peça
# m = 0.12125 + 0.4707 + 0.9238
# dm = 0.01*10^(-5) + 0.01*10^(-5) +  0.01*10^(-5)
# dm = 0.03*10^(-5)

global m = 1.51575;

global dm = 1e-5

global h = 0.467;
global dh = 1e-3;

global r = 6e-3;
global dr = 1e-4;

global g = 9.81;

global dt = 1e-2;

#=========================#
#  Conservação de energia #
#=========================#

function I(t)
    return ((g*(t^2))/(2*h) - 1)*m*r^2
end

function dI(t)
    termo1 = (2*r*dr*m + dm*r^2)*((g*t^2)/(2*h))
    termo2 = (g/2)*(((2*dt*t*h)+(dh*t^2))/h^2)*m*r^2
    termo2 = (2*r*dr*m + dm*r^2)
    return termo1 + termo2 - termo2
end

#==========================#
#  Pela geometria do corpo #
#==========================#

function Ic()
    raio = 6e-3
    M = 0.12125
    dM = 1e-5
    draio = 1e-4
    Ic = 0.5*(M)*(raio^2)
    dIc = raio*draio*M + 0.5*dM*(raio^2)
    println("Ic = " , Ic )
    println("ΔIc = ", dIc , "\n")
    return (Ic,dIc)
end

# id = Id()
# id[1] = valor da medida
# id[2] = valor da incerteza
function Id()
    raioInterno = 6e-4
    raioExterno = 6.25e-2

    M  = 0.4707
    dM = 1e-5
    draio = 1e-4

    Id = 0.5*M*(raioInterno^2 + raioExterno^2)
    dId = 0.5*dM*(raioInterno^2 + raioExterno^2) + M*(raioInterno*draio + raioExterno*draio)

 #  println("(1/2)*Δm*(r1^2 + r2^2) = ", dM*(raioInterno^2 + raioExterno^2))
 #  println("(1/2)*2*m*(r1*Δr + r2*Δr) = ", 0.5*2*M*(raioInterno*draio + raioExterno*draio))

    println("Id = ", Id)
    println("Δ Id = " , dId , "\n" )

    return (Id,dId)
end

function Ia()
    raioInterno = 0.0625
    raioExterno = 0.0760

    M  = 0.9238

    dM = 1e-5
    draio = 1e-4

    Ia = 0.5*M*(raioInterno^2 + raioExterno^2)
    dIa = 0.5*dM*(raioInterno^2 + raioExterno^2) + M*(raioInterno*draio + raioExterno*draio)

    println("Ia = ", Ia )
    println("Δ Ia = " , dIa , "\n"  )

    return (Ia,dIa)
end

function I()
    ia = Ia()
    ic = Ic()
    id = Id()
    println("I = Ia + Ic + Id \n I = " , ia[1] + ic[1] + id[1])
    println("\n ΔI = Δ Ia + ΔIc + ΔId \n Δ I= " , ia[2] + ic[2] + id[2])
end

max_stars_repo_path: Fisica experimental/02/01/src/exp1.jl

max_issues_repo_path: Fisica experimental/02/01/src/exp1.jl

max_forks_repo_path: Fisica experimental/02/01/src/exp1.jl


text:
const __bodyfunction__ = Dict{Method,Any}()

# Find keyword "body functions" (the function that contains the body
# as written by the developer, called after all missing keyword-arguments
# have been assigned values), in a manner that doesn't depend on
# gensymmed names.
# `mnokw` is the method that gets called when you invoke it without
# supplying any keywords.
function __lookup_kwbody__(mnokw::Method)
    function getsym(arg)
        isa(arg, Symbol) && return arg
        @assert isa(arg, GlobalRef)
        return arg.name
    end

    f = get(__bodyfunction__, mnokw, nothing)
    if f === nothing
        fmod = mnokw.module
        # The lowered code for `mnokw` should look like
        #   %1 = mkw(kwvalues..., #self#, args...)
        #        return %1
        # where `mkw` is the name of the "active" keyword body-function.
        ast = Base.uncompressed_ast(mnokw)
        if isa(ast, Core.CodeInfo) && length(ast.code) >= 2
            callexpr = ast.code[end-1]
            if isa(callexpr, Expr) && callexpr.head == :call
                fsym = callexpr.args[1]
                if isa(fsym, Symbol)
                    f = getfield(fmod, fsym)
                elseif isa(fsym, GlobalRef)
                    if fsym.mod === Core && fsym.name === :_apply
                        f = getfield(mnokw.module, getsym(callexpr.args[2]))
                    elseif fsym.mod === Core && fsym.name === :_apply_iterate
                        f = getfield(mnokw.module, getsym(callexpr.args[3]))
                    else
                        f = getfield(fsym.mod, fsym.name)
                    end
                else
                    f = missing
                end
            else
                f = missing
            end
        else
            f = missing
        end
        __bodyfunction__[mnokw] = f
    end
    return f
end

function _precompile_()
    ccall(:jl_generating_output, Cint, ()) == 1 || return nothing
    isdefined(Base.Broadcast, Symbol("#31#32")) && precompile(Tuple{getfield(Base.Broadcast, Symbol("#31#32")),DirectorySource})
    isdefined(Base.Broadcast, Symbol("#31#32")) && precompile(Tuple{getfield(Base.Broadcast, Symbol("#31#32")),FileSource})
    isdefined(BinaryBuilder, Symbol("#155#161")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#155#161")),String})
    isdefined(BinaryBuilder, Symbol("#170#173")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#170#173")),String})
    isdefined(BinaryBuilder, Symbol("#170#173")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#170#173")),String})
    isdefined(BinaryBuilder, Symbol("#20#21")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#20#21"))})
    isdefined(BinaryBuilder, Symbol("#20#21")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#20#21"))})
    isdefined(BinaryBuilder, Symbol("#222#227")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#222#227")),String})
    isdefined(BinaryBuilder, Symbol("#223#228")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#223#228")),String})
    isdefined(BinaryBuilder, Symbol("#224#229")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#224#229")),String})
    isdefined(BinaryBuilder, Symbol("#282#283")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#282#283"))})
    isdefined(BinaryBuilder, Symbol("#7#9")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#7#9"))})
    isdefined(BinaryBuilder, Symbol("#8#10")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#8#10"))})
    isdefined(BinaryBuilder, Symbol("#98#111")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#98#111")),String})
    isdefined(BinaryBuilder, Symbol("#99#112")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#99#112")),String})
    isdefined(BinaryBuilder, Symbol("#99#112")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#99#112")),SubString{String}})
    isdefined(BinaryBuilder, Symbol("#GOOS#208")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#GOOS#208")),Linux})
    isdefined(BinaryBuilder, Symbol("#base_gcc_flags#191")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#base_gcc_flags#191")),Linux})
    isdefined(BinaryBuilder, Symbol("#cargo#215")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#cargo#215")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#cc#205")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#cc#205")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#check_set#150")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#check_set#150")),String,Array{String,1}})
    isdefined(BinaryBuilder, Symbol("#clang#202")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#clang#202")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#clang_compile_flags#196")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#clang_compile_flags#196")),Linux})
    isdefined(BinaryBuilder, Symbol("#clang_flags#195")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#clang_flags#195")),Linux})
    isdefined(BinaryBuilder, Symbol("#clang_link_flags#197")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#clang_link_flags#197")),Linux})
    isdefined(BinaryBuilder, Symbol("#clangxx#203")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#clangxx#203")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#cxx#206")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#cxx#206")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#cxxfilt#217")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#cxxfilt#217")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#fortran_flags#193")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#fortran_flags#193")),Linux})
    isdefined(BinaryBuilder, Symbol("#gcc_flags#192")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#gcc_flags#192")),Linux})
    isdefined(BinaryBuilder, Symbol("#gcc_link_flags#198")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#gcc_link_flags#198")),Linux})
    isdefined(BinaryBuilder, Symbol("#gfortran#201")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#gfortran#201")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#go#210")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#go#210")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#meson#216")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#meson#216")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#objc#204")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#objc#204")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#rust_flags#212")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#rust_flags#212")),Linux})
    isdefined(BinaryBuilder, Symbol("#rustc#213")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#rustc#213")),IOStream,Linux})
    isdefined(BinaryBuilder, Symbol("#rustup#214")) && precompile(Tuple{getfield(BinaryBuilder, Symbol("#rustup#214")),IOStream,Linux})
    let fbody = try __lookup_kwbody__(which(BinaryBuilder.audit, (Prefix,String,))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Base.TTY,Linux,Bool,Bool,Bool,Bool,typeof(audit),Prefix,String,))
        end
    end
    let fbody = try __lookup_kwbody__(which(BinaryBuilder.build_tarballs, (Any,Any,Any,Any,Any,Any,Any,Any,))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Any,typeof(build_tarballs),Any,Any,Any,Any,Any,Any,Any,Any,))
        end
    end
    let fbody = try __lookup_kwbody__(which(BinaryBuilder.check_cxxstring_abi, (ObjectFile.ELF.ELFHandle{IOStream},Linux,))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Base.TTY,Bool,typeof(BinaryBuilder.check_cxxstring_abi),ObjectFile.ELF.ELFHandle{IOStream},Linux,))
        end
    end
    let fbody = try __lookup_kwbody__(which(BinaryBuilder.resolve_jlls, (Array{Dependency,1},))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Pkg.Types.Context,Base.TTY,typeof(BinaryBuilder.resolve_jlls),Array{Dependency,1},))
        end
    end
    let fbody = try __lookup_kwbody__(which(OutputCollector, (Cmd,))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Bool,Bool,Base.TTY,Type{OutputCollector},Cmd,))
        end
    end
    let fbody = try __lookup_kwbody__(which(run, (BinaryBuilder.UserNSRunner,Cmd,IOStream,))) catch missing end
        if !ismissing(fbody)
            precompile(fbody, (Bool,Base.TTY,typeof(run),BinaryBuilder.UserNSRunner,Cmd,IOStream,))
        end
    end
    precompile(Tuple{Core.kwftype(typeof(Base.Broadcast.broadcasted_kwsyntax)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(Base.Broadcast.broadcasted_kwsyntax),Function,Array{BinaryBuilder.AbstractSource,1}})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.Type)),NamedTuple{(:cwd, :platform),Tuple{String,Linux}},Type{BinaryBuilder.UserNSRunner},String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.Type)),NamedTuple{(:cwd, :platform, :verbose),Tuple{String,Linux,Bool}},Type{BinaryBuilder.UserNSRunner},String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.Type)),NamedTuple{(:cwd, :platform, :verbose, :workspaces, :compiler_wrapper_dir, :src_name, :shards, :compilers),Tuple{String,Linux,Bool,Array{Pair{String,String},1},String,String,Array{BinaryBuilder.CompilerShard,1},Array{Symbol,1}}},Type{BinaryBuilder.UserNSRunner},String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.autobuild)),NamedTuple{(:verbose, :debug, :meta_json_stream),Tuple{Bool,Bool,IOStream}},typeof(autobuild),String,String,VersionNumber,Array{DirectorySource,1},String,Array{Linux,1},Array{Product,1},Array{Dependency,1}})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.autobuild)),NamedTuple{(:verbose, :debug, :meta_json_stream, :compilers),Tuple{Bool,Bool,Nothing,Array{Symbol,1}}},typeof(autobuild),String,String,VersionNumber,Array{BinaryBuilder.AbstractSource,1},String,Array{Linux,1},Array{Product,1},Array{Dependency,1}})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.build_tarballs)),NamedTuple{(:compilers,),Tuple{Array{Symbol,1}}},typeof(build_tarballs),Array{String,1},String,VersionNumber,Array{BinaryBuilder.AbstractSource,1},String,Array{Linux,1},Array{Product,1},Array{Dependency,1}})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_cxxstring_abi)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(BinaryBuilder.check_cxxstring_abi),ObjectFile.ELF.ELFHandle{IOStream},Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_dynamic_linkage)),NamedTuple{(:platform, :silent, :verbose, :autofix),Tuple{Linux,Bool,Bool,Bool}},typeof(BinaryBuilder.check_dynamic_linkage),ObjectFile.ELF.ELFHandle{IOStream},Prefix,Array{String,1}})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_isa)),NamedTuple{(:verbose, :silent),Tuple{Bool,Bool}},typeof(BinaryBuilder.check_isa),ObjectFile.ELF.ELFHandle{IOStream},Linux,Prefix})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_libgfortran_version)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(BinaryBuilder.check_libgfortran_version),ObjectFile.ELF.ELFHandle{IOStream},Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_license)),NamedTuple{(:verbose, :silent),Tuple{Bool,Bool}},typeof(BinaryBuilder.check_license),Prefix,String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.check_os_abi)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(BinaryBuilder.check_os_abi),ObjectFile.ELF.ELFHandle{IOStream},Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.choose_shards)),NamedTuple{(:compilers,),Tuple{Array{Symbol,1}}},typeof(BinaryBuilder.choose_shards),Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.ensure_soname)),NamedTuple{(:verbose, :autofix),Tuple{Bool,Bool}},typeof(BinaryBuilder.ensure_soname),Prefix,String,Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.generate_compiler_wrappers!)),NamedTuple{(:bin_path, :compilers),Tuple{String,Array{Symbol,1}}},typeof(BinaryBuilder.generate_compiler_wrappers!),Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.generate_compiler_wrappers!)),NamedTuple{(:bin_path,),Tuple{String}},typeof(BinaryBuilder.generate_compiler_wrappers!),Linux})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.get_compilers_versions)),NamedTuple{(:compilers,),Tuple{Array{Symbol,1}}},typeof(BinaryBuilder.get_compilers_versions)})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.satisfied)),NamedTuple{(:verbose, :platform),Tuple{Bool,Linux}},typeof(satisfied),ExecutableProduct,Prefix})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.symlink_soname_lib)),NamedTuple{(:verbose, :autofix),Tuple{Bool,Bool}},typeof(BinaryBuilder.symlink_soname_lib),String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.translate_symlinks)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(BinaryBuilder.translate_symlinks),String})
    precompile(Tuple{Core.kwftype(typeof(BinaryBuilder.update_linkage)),NamedTuple{(:verbose,),Tuple{Bool}},typeof(BinaryBuilder.update_linkage),Prefix,Linux,String,SubString{String},String})
    precompile(Tuple{Type{Array{BinaryBuilder.AbstractSource,1}},UndefInitializer,Int64})
    precompile(Tuple{Type{Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Axes,F,Args} where Args<:Tuple where F where Axes},typeof(BinaryBuilder.coerce_dependency),Tuple{Array{Dependency,1}}})
    precompile(Tuple{Type{Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Axes,F,Args} where Args<:Tuple where F where Axes},typeof(BinaryBuilder.coerce_source),Tuple{Array{DirectorySource,1}}})
    precompile(Tuple{Type{Dependency},String})
    precompile(Tuple{Type{DirectorySource},String})
    precompile(Tuple{Type{ExecutableProduct},String,Symbol})
    precompile(Tuple{Type{FileSource},String,String})
    precompile(Tuple{Type{OutputCollector},Cmd,Base.Process,BinaryBuilder.LineStream,BinaryBuilder.LineStream,Base.GenericCondition{Base.AlwaysLockedST},Base.DevNull,Bool,Bool})
    precompile(Tuple{Type{OutputCollector},Cmd,Base.Process,BinaryBuilder.LineStream,BinaryBuilder.LineStream,Base.GenericCondition{Base.AlwaysLockedST},Base.TTY,Bool,Bool})
    precompile(Tuple{typeof(Base.Broadcast.broadcasted),Function,Array{BinaryBuilder.AbstractSource,1}})
    precompile(Tuple{typeof(Base.Broadcast.broadcasted),Function,Array{Dependency,1}})
    precompile(Tuple{typeof(Base.Broadcast.broadcasted),Function,Array{DirectorySource,1}})
    precompile(Tuple{typeof(Base.Broadcast.copyto_nonleaf!),Array{FileSource,1},Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(BinaryBuilder.coerce_source),Tuple{Base.Broadcast.Extruded{Array{BinaryBuilder.AbstractSource,1},Tuple{Bool},Tuple{Int64}}}},Base.OneTo{Int64},Int64,Int64})
    precompile(Tuple{typeof(Base.Broadcast.instantiate),Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Nothing,typeof(BinaryBuilder.coerce_dependency),Tuple{Array{Dependency,1}}}})
    precompile(Tuple{typeof(Base.Broadcast.instantiate),Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Nothing,typeof(BinaryBuilder.coerce_source),Tuple{Array{DirectorySource,1}}}})
    precompile(Tuple{typeof(Base.Broadcast.restart_copyto_nonleaf!),Array{BinaryBuilder.AbstractSource,1},Array{FileSource,1},Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(BinaryBuilder.coerce_source),Tuple{Base.Broadcast.Extruded{Array{BinaryBuilder.AbstractSource,1},Tuple{Bool},Tuple{Int64}}}},DirectorySource,Int64,Base.OneTo{Int64},Int64,Int64})
    precompile(Tuple{typeof(Base._compute_eltype),Type{Tuple{Pair{String,String},Pair{String,String},Pair{String,Array{DirectorySource,1}},Pair{String,String},Pair{String,Array{String,1}},Pair{String,Array{Product,1}},Pair{String,Array{Dependency,1}},Pair{String,Bool}}}})
    precompile(Tuple{typeof(Base.allocatedinline),Type{BinaryBuilder.CompilerShard}})
    precompile(Tuple{typeof(Base.grow_to!),Dict{String,Any},Tuple{Pair{String,String},Pair{String,String},Pair{String,Array{DirectorySource,1}},Pair{String,String},Pair{String,Array{String,1}},Pair{String,Array{Product,1}},Pair{String,Array{Dependency,1}},Pair{String,Bool}},Int64})
    precompile(Tuple{typeof(Base.grow_to!),Dict{String,String},Tuple{Pair{String,String},Pair{String,String},Pair{String,Array{DirectorySource,1}},Pair{String,String},Pair{String,Array{String,1}},Pair{String,Array{Product,1}},Pair{String,Array{Dependency,1}},Pair{String,Bool}},Int64})
    precompile(Tuple{typeof(Base.merge_types),NTuple{8,Symbol},Type{NamedTuple{(:cwd, :platform, :verbose, :workspaces, :compiler_wrapper_dir, :src_name, :shards),Tuple{String,Linux,Bool,Array{Pair{String,String},1},String,String,Array{BinaryBuilder.CompilerShard,1}}}},Type{NamedTuple{(:compilers,),Tuple{Array{Symbol,1}}}}})
    precompile(Tuple{typeof(Base.promote_typeof),FileSource,DirectorySource})
    precompile(Tuple{typeof(Base.vect),FileSource,Vararg{Any,N} where N})
    precompile(Tuple{typeof(BinaryBuilder.cppfilt),Array{SubString{String},1},Linux})
    precompile(Tuple{typeof(BinaryBuilder.minimum_instruction_set),Dict{String,Int64},Bool})
    precompile(Tuple{typeof(BinaryBuilder.should_ignore_lib),SubString{String},ObjectFile.ELF.ELFHandle{IOStream}})
    precompile(Tuple{typeof(BinaryBuilder.storage_dir),String,Vararg{String,N} where N})
    precompile(Tuple{typeof(BinaryBuilder.with_logfile),Function,Prefix,String})
    precompile(Tuple{typeof(JSON.Writer.show_element),JSON.Writer.CompactContext{Base.GenericIOBuffer{Array{UInt8,1}}},JSON.Serializations.StandardSerialization,ExecutableProduct})
    precompile(Tuple{typeof(JSON.Writer.show_pair),JSON.Writer.CompactContext{Base.GenericIOBuffer{Array{UInt8,1}}},JSON.Serializations.StandardSerialization,String,Array{Dependency,1}})
    precompile(Tuple{typeof(JSON.Writer.show_pair),JSON.Writer.CompactContext{Base.GenericIOBuffer{Array{UInt8,1}}},JSON.Serializations.StandardSerialization,String,Array{DirectorySource,1}})
    precompile(Tuple{typeof(JSON.Writer.show_pair),JSON.Writer.CompactContext{Base.GenericIOBuffer{Array{UInt8,1}}},JSON.Serializations.StandardSerialization,String,Array{Product,1}})
    precompile(Tuple{typeof(build_tarballs),Any,Any,Any,Any,Any,Any,Any,Any})
    precompile(Tuple{typeof(copy),Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(BinaryBuilder.coerce_dependency),Tuple{Array{Dependency,1}}}})
    precompile(Tuple{typeof(copy),Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(BinaryBuilder.coerce_source),Tuple{Array{DirectorySource,1}}}})
    precompile(Tuple{typeof(copyto!),Array{BinaryBuilder.AbstractSource,1},Tuple{FileSource,DirectorySource}})
    precompile(Tuple{typeof(findfirst),Function,Array{BinaryBuilder.CompilerShard,1}})
    precompile(Tuple{typeof(getindex),Type{Dependency},Dependency})
    precompile(Tuple{typeof(getindex),Type{Dependency}})
    precompile(Tuple{typeof(getindex),Type{DirectorySource}})
    precompile(Tuple{typeof(getindex),Type{Product},ExecutableProduct,ExecutableProduct,ExecutableProduct,ExecutableProduct,Vararg{ExecutableProduct,N} where N})
    precompile(Tuple{typeof(merge),NamedTuple{(:cwd, :platform, :verbose, :workspaces, :compiler_wrapper_dir, :src_name, :shards),Tuple{String,Linux,Bool,Array{Pair{String,String},1},String,String,Array{BinaryBuilder.CompilerShard,1}}},NamedTuple{(:compilers,),Tuple{Array{Symbol,1}}}})
    precompile(Tuple{typeof(repr),ExecutableProduct})
    precompile(Tuple{typeof(repr),LibraryProduct})
    precompile(Tuple{typeof(setindex!),Array{FileSource,1},FileSource,Int64})
    precompile(Tuple{typeof(setindex!),Dict{Any,Any},Dict{String,String},ExecutableProduct})
    precompile(Tuple{typeof(setindex!),Dict{String,Any},Array{Dependency,1},String})
    precompile(Tuple{typeof(setindex!),Dict{String,Any},Array{DirectorySource,1},String})
    precompile(Tuple{typeof(setindex!),Dict{String,Any},Array{Product,1},String})
    precompile(Tuple{typeof(similar),Base.Broadcast.Broadcasted{Base.Broadcast.DefaultArrayStyle{1},Tuple{Base.OneTo{Int64}},typeof(BinaryBuilder.coerce_source),Tuple{Base.Broadcast.Extruded{Array{BinaryBuilder.AbstractSource,1},Tuple{Bool},Tuple{Int64}}}},Type{FileSource}})
end


text:
module Utils
# utility functions

end #module


text:
using Metrics

@testset "Ranking and Statistical" begin

    y_rec = [0.85, 0.44, 0.71, 0.23, 0.90, 0.63, 0.16, 0.27]
    y_rel = [1, 0, 0, 1, 1, 1, 0, 0]
    output_rank_stats = Dict("precision_k" => 0.5, "recall_k" => 0.75, "f1_k" => 0.6)

    @testset "rank_stats_k" begin
        rank_out = ranking_stats_k(y_rec, y_rel; k=6)
        for i in keys(rank_out)
            @test rank_out[i] == output_rank_stats[i]
        end
    end

    @test avg_precision(y_rec, y_rel; k=6) == 0.6875
end


text:
module ImageGeoms

const RealU = Number # Union{Real, Unitful.Length}

abstract type Mask end
struct MaskAll <: Mask end
struct MaskCircle <: Mask end
struct MaskEllipse <: Mask end
struct MaskAllButEdge <: Mask end

export Mask, MaskAll, MaskCircle, MaskEllipse, MaskAllButEdge

include("core.jl")
include("imresize.jl")
include("makemask.jl")
include("mask.jl")

#include("property.jl") # exclude due to type inference issues

end # module


text:
using Test

include( "../src/Incomes.jl")

# using Incomes

@testset "Incomes" begin

    print( Incomes.isettostr(Incomes.EXEMPT_INCOME ))

end

max_stars_repo_name: grahamstark/ScottishTaxBenefitModel

max_issues_repo_name: grahamstark/ScottishTaxBenefitModel

max_forks_repo_name: grahamstark/ScottishTaxBenefitModel


text:
function translateToMathematica(x)
  rhs = """Set[{A021,A031,A032,A041,A042,A043,A121,A131,A132,A141,A142,A143,B021,B031,B032,B041,B042,B043,B121,B131,B132,B141,B142,B143,\\[Alpha]1,\\[Alpha]2,\\[Alpha]3,\\[Alpha]4,\\[Beta]11,\\[Beta]12,\\[Beta]13,\\[Beta]14,\\[Beta]21,\\[Beta]22,\\[Beta]23,\\[Beta]24,\\[Beta]31,\\[Beta]32,\\[Beta]33,\\[Beta]34,\\[Beta]41,\\[Beta]42,\\[Beta]43,\\[Beta]44}"""

  lhs = @sprintf("{%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s,%s}]",x[1],x[2],x[3],x[4],x[5],x[6],x[7],x[8],x[9],x[10],x[11],x[12],x[13],x[14],x[15],x[16],x[17],x[18],x[19],x[20],x[21],x[22],x[23],x[24],x[25],x[26],x[27],x[28],x[29],x[30],x[31],x[32],x[33],x[34],x[35],x[36],x[37],x[38],x[39],x[40],x[41],x[42],x[43],x[44])
  println(lhs)
  mathstr = "$rhs,$lhs\n"
  return mathstr
end

function MathematicaInitSet(x)
  lhs = @sprintf("{A021,%s},{A031,%s},{A032,%s},{A041,%s},{A042,%s},{A043,%s},{A121,%s},{A131,%s},{A132,%s},{A141,%s},{A142,%s},{A143,%s},{B021,%s},{B031,%s},{B032,%s},{B041,%s},{B042,%s},{B043,%s},{B121,%s},{B131,%s},{B132,%s},{B141,%s},{B142,%s},{B143,%s},{\\[Alpha]1,%s},{\\[Alpha]2,%s},{\\[Alpha]3,%s},{\\[Alpha]4,%s},{\\[Beta]11,%s},{\\[Beta]12,%s},{\\[Beta]13,%s},{\\[Beta]14,%s},{\\[Beta]21,%s},{\\[Beta]22,%s},{\\[Beta]23,%s},{\\[Beta]24,%s},{\\[Beta]31,%s},{\\[Beta]32,%s},{\\[Beta]33,%s},{\\[Beta]34,%s},{\\[Beta]41,%s},{\\[Beta]42,%s},{\\[Beta]43,%s},{\\[Beta]44,%s}\n",x[1],x[2],x[3],x[4],x[5],x[6],x[7],x[8],x[9],x[10],x[11],x[12],x[13],x[14],x[15],x[16],x[17],x[18],x[19],x[20],x[21],x[22],x[23],x[24],x[25],x[26],x[27],x[28],x[29],x[30],x[31],x[32],x[33],x[34],x[35],x[36],x[37],x[38],x[39],x[40],x[41],x[42],x[43],x[44])
  println(lhs)
  return lhs
end

function MathematicaReplaceSet(x)
  lhs = @sprintf("{A021->%s,A031->%s,A032->%s,A041->%s,A042->%s,A043->%s,A121->%s,A131->%s,A132->%s,A141->%s,A142->%s,A143->%s,B021->%s,B031->%s,B032->%s,B041->%s,B042->%s,B043->%s,B121->%s,B131->%s,B132->%s,B141->%s,B142->%s,B143->%s,\\[Alpha]1->%s,\\[Alpha]2->%s,\\[Alpha]3->%s,\\[Alpha]4->%s,\\[Beta]11->%s,\\[Beta]12->%s,\\[Beta]13->%s,\\[Beta]14->%s,\\[Beta]21->%s,\\[Beta]22->%s,\\[Beta]23->%s,\\[Beta]24->%s,\\[Beta]31->%s,\\[Beta]32->%s,\\[Beta]33->%s,\\[Beta]34->%s,\\[Beta]41->%s,\\[Beta]42->%s,\\[Beta]43->%s,\\[Beta]44->%s}\n",x[1],x[2],x[3],x[4],x[5],x[6],x[7],x[8],x[9],x[10],x[11],x[12],x[13],x[14],x[15],x[16],x[17],x[18],x[19],x[20],x[21],x[22],x[23],x[24],x[25],x[26],x[27],x[28],x[29],x[30],x[31],x[32],x[33],x[34],x[35],x[36],x[37],x[38],x[39],x[40],x[41],x[42],x[43],x[44])
  println(lhs)
  return lhs
end

function printForJulia(x)
  A0,A1,B0,B1,α,β1,β2,β3,β4 = translate(x)
  eV = ones(4)
  c0 = A0*ones(4)
  c1 = A1*ones(4)
  resString = """
  x = [$(join(string.(x), ","))]
  c0 = [$(join(string.(c0), ","))]
  c1 = [$(join(string.(c1), ","))]
  A0 = [$(join(string.(A0), ","))]
  A1 = [$(join(string.(A1), ","))]
  B0 = [$(join(string.(B0), ","))]
  B1 = [$(join(string.(B1), ","))]
  α  = [$(join(string.(α), ","))]
  β1 = [$(join(string.(β1), ","))]
  β2 = [$(join(string.(β2), ","))]
  β3 = [$(join(string.(β3), ","))]
  β4 = [$(join(string.(β4), ","))]
  """
  return resString
end


text:
#nb # %% slideshow [markdown] {"slideshow": {"slide_type": "slide"}}
#md # # [Lecture 5 - Container types and indexing](@id lecture5)
#!md # # Lecture 5 - Container types and indexing
#
# **Essential Skills For Computational Biology**
#
# *2020-06-16*
#
# Kevin Bonham, PhD

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "slide"}}
# ## Scalars and Containers

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# - A "scalar" is a singular value

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# - Most of the types you've seen so far are scalars
#   - Eg. `Float64`, `Int64`, `Bool`, `Char`

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# - There are also many container types, which store multiple values

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# - Containers can usually be accessed by an "index"

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## A Tuple is a simple, ordered container

my_tuple = (1, 1.2, "a string")

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# Each position has an integer index,
# 1 for the first position, 2 for the second, and so on
my_tuple[1]
#-
my_tuple[3]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Attempting to access an index that doesn't exist is an error
length(my_tuple)
#-
my_tuple[4]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "slide"}}
# ## A Vector is a mutable, ordered container
#
# Vectors can by easily created with `[` and `]`

#nb %% A slide [code] {"slideshow": {"slide_type": "fragment"}}
my_vector = [1,3,5,7]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# Each position has an integer index,
# 1 for the first position, 2 for the second, and so on

#nb %% A slide [code] {"slideshow": {"slide_type": "fragment"}}
my_vector[1]
#-
my_vector[2]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Add to the end of a `Vector` with `push!()`

#nb %% A slide [code] {"slideshow": {"slide_type": "fragment"}}
push!(my_vector, 50);
#-
my_vector

#nb %% A slide [code] {"slideshow": {"slide_type": "subslide"}}
other_vector = [] ## an empty vector
for i in 1:2:10
    push!(other_vector, i)
end
#-
other_vector

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Change the value of a vector by assignment (`=`)

#nb %% A slide [code] {"slideshow": {"slide_type": "fragment"}}
my_vector[1] = 42;
#-
my_vector

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Index a "slice" with a range
my_vector
#-
my_vector[2:4]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Attempting to access an index that doesn't exist is an error
length(my_vector)
#-
my_vector[6]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "slide"}}
# # Dictionaries are containers with "key" => "value" pairs
#
# - Keys can be any scalar type
# - Keys are used as indexes
# - Values can be scalars or other containers
# - Keys are unordered

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Create a Dictionary with `Dict()`
my_dict = Dict("key1"=> 42, "key2"=> 1.8, 5=> "a value")

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## The keys are used to index
my_dict["key1"]
#-
my_dict[5]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Values can be reassigned
my_dict["key1"] = 3.14

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "subslide"}}
# ## Attempting to access with a key that doesn't exist is an error
my_dict["I don't exist"]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# But you can assign a new key-value pair
my_dict["I don't exist... yet"] = "Now I do!"
#-
my_dict["I don't exist... yet"]

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "fragment"}}
# You can use `haskey()` (a Boolean function) to see if a dictionary already has a key
haskey(my_dict, "key1")
#-
haskey(my_dict, "key1000")

#nb # %% A slide [markdown] {"slideshow": {"slide_type": "slide"}}
# # Questions?