text:
###
# Authors: Daniel Thuerck
# NEC Laboratories Europe GmbH, Copyright (c) 2021, All rights reserved.
#
# THIS HEADER MAY NOT BE EXTRACTED OR MODIFIED IN ANY WAY.
#
# SOFTWARE LICENSE AGREEMENT
# ACADEMIC OR NON-PROFIT ORGANIZATION NONCOMMERCIAL RESEARCH USE ONLY
# BY USING OR DOWNLOADING THE SOFTWARE, YOU ARE AGREEING TO THE TERMS
# OF THIS LICENSE AGREEMENT. IF YOU DO NOT AGREE WITH THESE TERMS, YOU MAY
# NOT USE OR DOWNLOAD THE SOFTWARE.
#
# This is a license agreement ("Agreement") between your academic institution or non-profit
# organization or self (called "Licensee" or "You" in this Agreement) and NEC Laboratories Europe
# GmbH (called "Licensor" in this Agreement). All rights not specifically granted to you in this
# Agreement are reserved for Licensor.
#
# RESERVATION OF OWNERSHIP AND GRANT OF LICENSE: Licensor retains exclusive
# ownership of any copy of the Software (as defined below) licensed under this Agreement and hereby
# grants to Licensee a personal, non-exclusive, non-transferable license to use the Software for
# noncommercial research purposes, without the right to sublicense, pursuant to the terms and
# conditions of this Agreement. NO EXPRESS OR IMPLIED LICENSES TO ANY OF LICENSOR’S
# PATENT RIGHTS ARE GRANTED BY THIS LICENSE. As used in this Agreement, the term
# "Software" means (i) the actual copy of all or any portion of code for program routines made
# accessible to Licensee by Licensor pursuant to this Agreement, inclusive of backups, updates, and/or
# merged copies permitted hereunder or subsequently supplied by Licensor, including all or any file
# structures, programming instructions, user interfaces and screen formats and sequences as well as any
# and all documentation and instructions related to it, and (ii) all or any derivatives and/or modifications
# created or made by You to any of the items specified in (i).
#
# CONFIDENTIALITY/PUBLICATIONS: Licensee acknowledges that the Software is proprietary to
# Licensor, and as such, Licensee agrees to receive all such materials and to use the Software only in
# accordance with the terms of this Agreement. Licensee agrees to use reasonable effort to protect the
# Software from unauthorized use, reproduction, distribution, or publication. All publication materials
# mentioning features or use of this software must explicitly include an acknowledgement the software
# was developed by NEC Laboratories Europe GmbH.
#
# COPYRIGHT: The Software is owned by Licensor.
#
# PERMITTED USES: The Software may be used for your own noncommercial internal research
# purposes. You understand and agree that Licensor is not obligated to implement any suggestionsand/or
# feedback you might provide regarding the Software, but to the extent Licensor does so, you are
# not entitled to any compensation related thereto.
#
# DERIVATIVES: You may create derivatives of or make modifications to the Software, however, You
# agree that all and any such derivatives and modifications will be owned by Licensor and become a part
# of the Software licensed to You under this Agreement. You may only use such derivatives and
# modifications for your own noncommercial internal research purposes, and you may not otherwise
# use, distribute or copy such derivatives and modifications in violation of this Agreement.
#
# BACKUPS: If Licensee is an organization, it may make that number of copies of the Software
# necessary for internal noncommercial use at a single site within its organization provided that all
# information appearing in or on the original labels, including the copyright and trademark notices are
# copied onto the labels of the copies.
#
# USES NOT PERMITTED: You may not distribute, copy or use the Software except as explicitly
# permitted herein. Licensee has not been granted any trademark license as part of this Agreement.
# Neither the name of NEC Laboratories Europe GmbH nor the names of its contributors may be used to
# endorse or promote products derived from this Software without specific prior written permission.
# You may not sell, rent, lease, sublicense, lend, time-share or transfer, in whole or in part, or provide
# third parties access to prior or present versions (or any parts thereof) of the Software.
#
# ASSIGNMENT: You may not assign this Agreement or your rights hereunder without the prior
# written consent of Licensor. Any attempted assignment without such consent shall be null and void.
#
# TERM: The term of the license granted by this Agreement is from Licensee's acceptance of this
# Agreement by downloading the Software or by using the Software until terminated as provided below.
# The Agreement automatically terminates without notice if you fail to comply with any provision of
# this Agreement. Licensee may terminate this Agreement by ceasing using the Software. Upon any
# termination of this Agreement, Licensee will delete any and all copies of the Software. You agree that
# all provisions which operate to protect the proprietary rights of Licensor shall remain in force should
# breach occur and that the obligation of confidentiality described in this Agreement is binding in
# perpetuity and, as such, survives the term of the Agreement.
#
# FEE: Provided Licensee abides completely by the terms and conditions of this Agreement, there is no
# fee due to Licensor for Licensee's use of the Software in accordance with this Agreement.
#
# DISCLAIMER OF WARRANTIES: THE SOFTWARE IS PROVIDED "AS-IS" WITHOUT
# WARRANTY OF ANY KIND INCLUDING ANY WARRANTIES OF PERFORMANCE OR
# MERCHANTABILITY OR FITNESS FOR A PARTICULAR USE OR PURPOSE OR OF NON-
# INFRINGEMENT. LICENSEE BEARS ALL RISK RELATING TO QUALITY AND
# PERFORMANCE OF THE SOFTWARE AND RELATED MATERIALS.
#
# SUPPORT AND MAINTENANCE: No Software support or training by the Licensor is provided as
# part of this Agreement.
#
# EXCLUSIVE REMEDY AND LIMITATION OF LIABILITY: To the maximum extent permitted
# under applicable law, Licensor shall not be liable for direct, indirect, special, incidental, or
# consequential damages or lost profits related to Licensee's use of and/or inability to use the Software,
# even if Licensor is advised of the possibility of such damage.
#
# EXPORT REGULATION: Licensee agrees to comply with any and all applicable export control laws,
# regulations, and/or other laws related to embargoes and sanction programs administered by law.
#
# SEVERABILITY: If any provision(s) of this Agreement shall be held to be invalid, illegal, or
# unenforceable by a court or other tribunal of competent jurisdiction, the validity, legality and
# enforceability of the remaining provisions shall not in any way be affected or impaired thereby.
#
# NO IMPLIED WAIVERS: No failure or delay by Licensor in enforcing any right or remedy under this
# Agreement shall be construed as a waiver of any future or other exercise of such right or remedy by
# Licensor.
#
# GOVERNING LAW: This Agreement shall be construed and enforced in accordance with the laws of
# Germany without reference to conflict of laws principles. You consent to the personal jurisdiction of
# the courts of this country and waive their rights to venue outside of Germany.
#
# ENTIRE AGREEMENT AND AMENDMENTS: This Agreement constitutes the sole and entire
# agreement between Licensee and Licensor as to the matter set forth herein and supersedes any
# previous agreements, understandings, and arrangements between the parties relating hereto.
#
# THIS HEADER MAY NOT BE EXTRACTED OR MODIFIED IN ANY WAY.
###

include("typestructs.jl")

###
# General tools
###
primitive_cmap = Dict(
    :Nothing => Cvoid,
    :Number => Cdouble,
	:Complex => ComplexF64,
	:Real => Cdouble,
    :AbstractFloat => Cdouble,
    :BigFloat => Cdouble,
    :Float16 => Cfloat,
    :Float32 => Cfloat,
    :Float64 => Cdouble,
    :Integer => Clonglong,
    :BigInt => Clonglong,
    :Bool => Cuchar,
    :Signed => Clonglong,
    :Int128 => Clonglong,
    :Int16 => Cshort,
    :Int32 => Cint,
    :Int64 => Clonglong,
    :Int8 => Cchar,
    :Unsigned => Culonglong,
    :UInt128 => Culonglong,
    :UInt16 => Cushort,
    :UInt32 => Cuint,
    :UInt64 => Culonglong,
    :UInt8 => Cuchar
)

# creating unmanaged type (e.g. ptr + dims for arrays)
function create_unmanaged_type_replacements(t::param, reg :: type_definition_registry)
    t.managed = false

    if isa(t.type, struct_ref_t)
        t.managed = true

        # just replace by pointer
        t.replacement = param(
            t.name,
            ptr_t(t.type),
            false,
            nothing)
    end

    if isa(t.type, array_t)
        t.managed = true

        # replace by pointer to an equivalent struct
        t_sym = Symbol(hash_type(t.type))
        if !(t_sym in keys(reg.defs))
            reg.defs[t_sym] = t.type
            reg.names[t_sym] = Symbol("jl_implicit_" * string(length(keys(reg.defs))))
        end

        t_sym = reg.names[t_sym]

        tpl_struct = struct_ref_t(t_sym)

        t.replacement = param(
            t.name,
            ptr_t(tpl_struct),
            false,
            nothing)
    end

    if isa(t.type, tuple_t)
        t.managed = true

        # replace by pointer to an equivalent struct
        t_sym = Symbol(hash_type(t.type))
        if !(t_sym in keys(reg.defs))
            reg.defs[t_sym] = t.type
            reg.names[t_sym] = Symbol("jl_implicit_" * string(length(keys(reg.defs))))
        end

        t_sym = reg.names[t_sym]

        tpl_struct = struct_ref_t(t_sym)

        t.replacement = param(
            t.name,
            ptr_t(tpl_struct),
            false,
            nothing)
    end

end

function generate_implicit_structs(reg::type_definition_registry)

    for key in keys(reg.defs)
        stpl = reg.defs[key]
        members = []

        if isa(stpl, tuple_t)
            for ix in 1:length(stpl.members)
                push!(
                    members,
                    param(
                        Symbol("elem" * string(ix)),
                        stpl.members[ix],
                        false,
                        [])
                    )
            end
        elseif isa(stpl, array_t)
            # replace by pointer to native type and dimensions
            push!(
                members,
                param(
                    Symbol("ptr"),
                    ptr_t(stpl.pointee_t),
                    false,
                    []
                ))

            for d in 1 : stpl.dims
                push!(
                    members,
                    param(
                        Symbol("dim" * string(d)),
                        integer_t(64, true),
                        false,
                        []
                    ))
            end
        end

        push!(reg.adapted, (reg.names[key], members))
    end
end

###
# Julia
###

function deserialize(type_reg::type_definition_registry, val_type::primitive_t, replace_struct_t, replace_struct_ptr, val_expr)
    ins = [Expr(
        Symbol("="),
        val_expr,
        Expr(
        Symbol("call"),
        Symbol("unsafe_load"),
        Expr(
            Symbol("call"),
            Symbol("reinterpret"),
            expr_from_type(ptr_t(val_type)),
            replace_struct_ptr)))
    ]

    return ins
end

function deserialize(type_reg::type_definition_registry, val_type::ptr_t, replace_struct_t, replace_struct_ptr, val_expr)
    ins = [Expr(
        Symbol("="),
        val_expr,
        Expr(
        Symbol("call"),
        Symbol("unsafe_load"),
        Expr(
            Symbol("call"),
            Symbol("reinterpret"),
            expr_from_type(ptr_t(val_type)),
            replace_struct_ptr)))
    ]

    return ins
end

function deserialize(type_reg::type_definition_registry, val_type::array_t, replace_struct_t, replace_struct_ptr, val_expr)
    # array = structure with value pointer and one int per dimension
    ins = []

    # dimensions
    dim_expr = nothing

    dim_list = []
    for dim in 1:val_type.dims
        dim_sym = Symbol(string(val_expr) * "_dim" * string(dim))

        push!(
            ins,
            Expr(
                Symbol("="),
                dim_sym,
                generate_struct_field_read(
                    replace_struct_ptr,
                    replace_struct_t,
                    integer_t(64, false),
                    1 + dim
                )
            )
        )
        push!(dim_list, dim_sym)
    end

    if val_type.dims > 1
        dim_expr = Expr(
            Symbol("tuple"),
            dim_list...
        )
    else
        dim_expr = dim_list[1]
    end

    # create Array{...} from data -> stored in val_expr
    push!(
        ins,
        Expr(
            Symbol("="),
            val_expr,
            Expr(
                Symbol("call"),
                Symbol("unsafe_wrap"),
                expr_from_type(val_type),
                generate_struct_field_read(
                    replace_struct_ptr,
                    replace_struct_t,
                    ptr_t(val_type.pointee_t),
                    1
                ),
                dim_expr
            )
        )
    )

    # push!(
    #     ins,
    #     Expr(
    #         Symbol("="),
    #         val_expr,
    #         Expr(
    #             Symbol("call"),
    #             Symbol("unsafe_wrap"),
    #             expr_from_type(val_type),
    #             Expr(
    #                 Symbol("call"),
    #                 Symbol("unsafe_load"),
    #                 struct_member_ptr(
    #                     replace_struct_ptr,
    #                     replace_struct_t,
    #                     ptr_t(val_type.pointee_t),
    #                     1
    #                 )
    #             ),
    #             dim_expr
    #         )
    #     )
    # )

    return ins
end

function deserialize(type_reg::type_definition_registry, val_type::struct_ref_t, replace_struct_t, replace_struct_ptr, val_expr)
    ins = []
    ctr = 0

    # retrieve struct definition (note: struct_ref_t as original type always refers
    # to explicit structs -> hence, there is no replacement)
    # note that once struct members have been deserialized, we need to instantiate the
    # original structure, not the replacement that we have created
    hash_name = val_type.name.name
    val_type_def = type_reg.defs[val_type.name.name]

    # println("Deserialize struct: ", hash_name, " with replacement type ", replace_struct_t)

    # deserialize each individual member
    member_vars = []
    for (ix, member) in enumerate(val_type_def.members)
        member_var = Symbol(string(val_expr) * "_" * string(ctr))

        # println(ix, " -> ", member)

        # compute offset to this member
        member_ptr = nothing
        if isa(member.type, primitive_t)
            member_ptr = struct_member_ptr(
                replace_struct_ptr,
                replace_struct_t,
                member.type,
                ix)
        else
            member_ptr = generate_struct_field_read(
                replace_struct_ptr,
                replace_struct_t,
                ptr_t(member.replacement),
                ix)
        end

        append!(
            ins,
            deserialize(
                type_reg,
                member.type,
                isa(member.type, primitive_t) ? nothing : member.replacement,
                member_ptr,
                member_var
            )
        )

        push!(member_vars, member_var)

        ctr += 1
    end

    # now deserialize the _actual_ structure
    push!(
        ins,
        Expr(
            Symbol("="),
            val_expr,
            Expr(
                Symbol("call"),
                val_type_def.name.fqname,
                member_vars...
            )
        )
    )

    return ins
end

function deserialize(type_reg::type_definition_registry, val_type::tuple_t, replace_struct_t, replace_struct_ptr, val_expr)
    ins = []
    ctr = 0

    # retrieve struct definition
    val_type_def = type_reg.defs[replace_struct_t.name.name]

    # println("Deserialize tuple: ", val_type, " with replacement type ", replace_struct_t)

    # deserialize each individual member
    member_vars = []
    for (ix, member) in enumerate(val_type_def.members)
        member_var = Symbol(string(val_expr) * "_" * string(ctr))

        # println("-> ", member)

        # compute offset to this member
        member_ptr = nothing
        if isa(member.type, primitive_t)
            member_ptr = struct_member_ptr(
                replace_struct_ptr,
                replace_struct_t,
                member.type,
                ix,
                true)
        else
            member_ptr = generate_struct_field_read(
                replace_struct_ptr,
                replace_struct_t,
                ptr_t(member.replacement),
                ix)
        end

        append!(
            ins,
            deserialize(
                type_reg,
                member.type,
                isa(member.type, primitive_t) ? nothing : member.replacement,
                member_ptr,
                member_var
            )
        )

        push!(member_vars, member_var)

        ctr += 1
    end

    # now deserialize the _actual_ structure
    push!(
        ins,
        Expr(
            Symbol("="),
            val_expr,
            Expr(
                Symbol("tuple"),
                member_vars...
            )
        )
    )

    return ins
end

function serialize!(ptr_sym, ret_sym, ret_type::array_t, repl_struct_ref, type_defs::type_definition_registry)
    ins = []

    # alloc pointer
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Expr(
                    Symbol("."),
                    Symbol("Libc"),
                    QuoteNode(Symbol("malloc"))
                ),
                Expr(
                    Symbol("call"),
                    Symbol("sizeof"),
                    repl_struct_ref.name.name
                )
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Symbol("reinterpret"),
                Expr(
                    Symbol("curly"),
                    Symbol("Ptr"),
                    repl_struct_ref.name.name
                ),
                ptr_sym
            )
        )
    )

    # compute number of elements
    dim_exprs = []
    for dim in 1:ret_type.dims
        push!(
            dim_exprs,
            Expr(
                Symbol("call"),
                Symbol("size"),
                ret_sym,
                dim
            )
        )
    end

    dim_expr = Expr(
        Symbol("call"),
        Symbol("size"),
        ret_sym,
        1
    )

    for dim in 2:ret_type.dims
        dim_expr = Expr(
            Symbol("call"),
            Symbol("*"),
            dim_expr,
            Expr(
                Symbol("call"),
                Symbol("size"),
                ret_sym,
                dim
            )
        )
    end

    # write data to storage (after allocating memory for it)
    val_ptr_sym = Symbol(String(ptr_sym) * "_arr")
    push!(
        ins,
        Expr(
            Symbol("="),
            val_ptr_sym,
            Expr(
                Symbol("call"),
                Expr(
                    Symbol("."),
                    Symbol("Libc"),
                    QuoteNode(Symbol("malloc"))
                ),
                Expr(
                    Symbol("call"),
                    Symbol("*"),
                    dim_expr,
                    Expr(
                        Symbol("call"),
                        Symbol("sizeof"),
                        expr_from_type(ret_type.pointee_t)
                    )
                )
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("="),
            val_ptr_sym,
            Expr(
                Symbol("call"),
                Symbol("reinterpret"),
                Expr(
                    Symbol("curly"),
                    Symbol("Ptr"),
                    expr_from_type(ret_type.pointee_t)
                ),
                val_ptr_sym
            )
        )
    )

    push!(
        ins,
        Expr(
            Symbol("call"),
            Symbol("unsafe_copyto!"),
            val_ptr_sym,
            Expr(
                Symbol("call"),
                Symbol("pointer"),
                ret_sym
            ),
            dim_expr
        )
    )

    # create return struct,, then serialize it
    ret_val = Symbol(String(ptr_sym) * "_val")
    push!(
        ins,
        Expr(
            Symbol("="),
            ret_val,
            Expr(
                Symbol("call"),
                repl_struct_ref.name.name,
                val_ptr_sym,
                dim_exprs...
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("call"),
            Symbol("unsafe_store!"),
            ptr_sym,
            ret_val
        )
    )

    return ins
end

function serialize!(ptr_sym, ret_sym, ret_type::struct_ref_t, repl_struct_ref, type_defs::type_definition_registry)
    ins = []

    # alloc pointer
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Expr(
                    Symbol("."),
                    Symbol("Libc"),
                    QuoteNode(Symbol("malloc"))
                ),
                Expr(
                    Symbol("call"),
                    Symbol("sizeof"),
                    repl_struct_ref.name.name
                )
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Symbol("reinterpret"),
                Expr(
                    Symbol("curly"),
                    Symbol("Ptr"),
                    repl_struct_ref.name.name
                ),
                ptr_sym
            )
        )
    )

    # serialize members
    member_exprs = []

    val_type_def = type_defs.defs[ret_type.name.name]

    for (ix, member) in enumerate(val_type_def.members)
        acc_expr = Expr(
            Symbol("."),
            ret_sym,
            QuoteNode(member.name)
        )

        if isa(member.type, primitive_t)
            # primitive value -> just copy
            push!(member_exprs, acc_expr)
        else
            # complex type -> deserialize to additional pointer
            member_sym = Symbol(String(ptr_sym) * "_" * string(ix))
            append!(
                ins,
                serialize!(member_sym, acc_expr, member.type, member.replacement, type_defs)
            )
            push!(member_exprs, member_sym)
        end
    end

    # create return struct,, then serialize it
    ret_val = Symbol(String(ptr_sym) * "_val")
    push!(
        ins,
        Expr(
            Symbol("="),
            ret_val,
            Expr(
                Symbol("call"),
                repl_struct_ref.name.name,
                member_exprs...
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("call"),
            Symbol("unsafe_store!"),
            ptr_sym,
            ret_val
        )
    )

end

function serialize!(ptr_sym, ret_sym, ret_type::tuple_t, repl_struct_ref, type_defs::type_definition_registry)
    ins = []

    # alloc pointer
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Expr(
                    Symbol("."),
                    Symbol("Libc"),
                    QuoteNode(Symbol("malloc"))
                ),
                Expr(
                    Symbol("call"),
                    Symbol("sizeof"),
                    repl_struct_ref.name.name
                )
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("="),
            ptr_sym,
            Expr(
                Symbol("call"),
                Symbol("reinterpret"),
                Expr(
                    Symbol("curly"),
                    Symbol("Ptr"),
                    repl_struct_ref.name.name
                ),
                ptr_sym
            )
        )
    )

    # serialize members
    member_exprs = []

    val_type_def = type_defs.defs[repl_struct_ref.name.name]

    for (ix, member) in enumerate(val_type_def.members)
        acc_expr = Expr(
            Symbol("ref"),
            ret_sym,
            ix
        )

        if isa(member.type, primitive_t)
            # primitive value -> just copy
            push!(member_exprs, acc_expr)
        else
            # complex type -> deserialize to additional pointer
            member_sym = Symbol(String(ptr_sym) * "_" * string(ix))
            append!(
                ins,
                serialize!(member_sym, acc_expr, member.type, member.replacement, type_defs)
            )
            push!(member_exprs, member_sym)
        end
    end

    # create return struct,, then serialize it
    ret_val = Symbol(String(ptr_sym) * "_val")
    push!(
        ins,
        Expr(
            Symbol("="),
            ret_val,
            Expr(
                Symbol("call"),
                repl_struct_ref.name.name,
                member_exprs...
            )
        )
    )
    push!(
        ins,
        Expr(
            Symbol("call"),
            Symbol("unsafe_store!"),
            ptr_sym,
            ret_val
        )
    )

end

function generate_julia_struct_def(name, members)
    # create list of members
    expr_members = []

    for member in members
        push!(
            expr_members,
            Expr(
                Symbol("::"),
                member.name,
                expr_from_type(member.managed ? ptr_t(member.replacement) : member.type)
            )
        )
    end

    return [
        Expr(
            Symbol("struct"),
            true,
            name,
            Expr(
                Symbol("block"),
                expr_members...
            )
        ),
        Expr(
            Symbol("export"),
            Symbol(name))
    ]
end

function generate_struct_field_write(struct_ptr, struct_ref_type, elem_type, member_ix, rhs_expr)
    return Expr(
        Symbol("call"),
        Symbol("unsafe_store!"),
        Expr(
            Symbol("call"),
            Symbol("reinterpret"),
            expr_from_type(ptr_t(elem_type)),
            Expr(
                Symbol("call"),
                Symbol("+"),
                struct_ptr,
                Expr(
                    Symbol("call"),
                    Symbol("fieldoffset"),
                    expr_from_type(struct_ref_type),
                    member_ix)
            )
        ),
        rhs_expr
    )
end

function struct_member_ptr(struct_ptr, struct_ref_type, elem_type, member_ix, raw = false)
    raw_ptr = Expr(
        Symbol("call"),
        Symbol("+"),
        Expr(
            Symbol("call"),
            Symbol("reinterpret"),
            expr_from_type(ptr_t(integer_t(8, false))),
            struct_ptr
        ),
        Expr(
            Symbol("call"),
            Symbol("fieldoffset"),
            expr_from_type(struct_ref_type),
            member_ix)
    )

    if raw
        return raw_ptr
    end

    return Expr(
        Symbol("call"),
        Symbol("reinterpret"),
        expr_from_type(ptr_t(elem_type)),
        raw_ptr
    )
end

function generate_struct_field_read(struct_ptr, struct_ref_type, elem_type, member_ix)
    return Expr(
        Symbol("call"),
        Symbol("unsafe_load"),
        struct_member_ptr(struct_ptr, struct_ref_type, elem_type, member_ix)
    )
end

function generate_julia_wrapper_function(type_reg::type_definition_registry, name, call_sym, ret, params)
    fun_name = Symbol("jl_" * string(name))

    # function name (plus prefix)
    call_args = []
    push!(call_args, fun_name)

    # arguments -> use pointers for structs
    for param in params
        # println(param)

        push!(
            call_args,
            Expr(
                Symbol("::"),
                param.name,
                expr_from_type(param.replacement === nothing ? param.type : ptr_t(param.replacement))
            )
        )
    end

    signature = Expr(
        Symbol("::"),
        Expr(
            :call,
            call_args...
        ),
        expr_from_type(ret.replacement === nothing ? ret.type : ptr_t(ret.replacement))
    )

    ##
    # deserialization
    ##
    deserialized_args = []
    ins = []

    for (ctr, param) in enumerate(params)
        # println("PARAM ", param.name, " ", param.type, " ", isa(param.type, ptr_t) ? param.type.pointee_t : "")

        if isa(param.type, primitive_t)
            push!(
                deserialized_args,
                param.name)
        else
            param_var = Symbol("var_" * string(ctr))

            append!(
                ins,
                deserialize(
                    type_reg,
                    param.type,
                    param.replacement,
                    param.name,
                    param_var
                )
            )
            push!(
                deserialized_args,
                param_var
            )
        end
    end

    ##
    # call actual function
    ##
    call_expr = Expr(
        Symbol("call"),
        # Symbol(name),
        call_sym,
        deserialized_args...
    )

    ret_val_sym = Symbol("ret_val")
    if typeof(ret.type) != void_t
        if isa(ret.type, primitive_t)
            call_expr = Expr(
                Symbol("return"),
                call_expr
            )
        else
            call_expr = Expr(
                Symbol("="),
                ret_val_sym,
                call_expr
            )
        end
    end

    push!(
        ins,
        call_expr
    )

    ##
    # serialization of return values
    ##
    wrap_in_return = false
    if typeof(ret.type) != void_t
        if !isa(ret.type, primitive_t)
            ret_ptr_sym = Symbol("ret_ptr")
            append!(
                ins,
                serialize!(ret_ptr_sym, ret_val_sym, ret.type, ret.replacement,
                    type_reg)
            )

            # complex return type, needs to be serialized to
            # memory that is allocated by Julia's Libc-Wrapper

            # return the pointer to C
            push!(
                ins,
                Expr(
                    Symbol("return"),
                    ret_ptr_sym
                ))
        end
    end

    ##
    # macrocall structure
    ##
    return [
        Expr(
            :macrocall,
            Expr(
                Symbol("."),
                :Base,
                QuoteNode(Symbol("@ccallable"))
            ),
            Nothing,
            Expr(
                :function,
                signature,
                Expr(
                    :block,
                    ins...
                )
            )
        ),
        Expr(
            Symbol("export"),
            fun_name
        )
    ]
end

###
# C
###

function generate_c_struct_def(name, members)
    cstr = "typedef struct {\n"

    # generate a list of unmanaged members
    unm_members = []
    for p in members
        if p.managed
            push!(unm_members, param(p.name, ptr_t(p.replacement), false, nothing))
        elseif typeof(p.type) == struct_ref_t
            push!(unm_members, param(p.name, ptr_t(p.type), false, nothing))
        else
            push!(unm_members, p)
        end
    end

    for member in unm_members
        cstr *= "\t" * cstring_from_type(member.type) * " " * string(member.name) * ";\n"
    end
    cstr *= "} " * string(name) * ";\n"

    return cstr
end

function generate_c_function_def(name, rettype, params)
    cstr = cstring_from_type(rettype.replacement === nothing ? rettype.type :
        ptr_t(rettype.replacement))
    cstr *= " jl_" * string(name) * "(\n"

    # generate a list of unmanaged members
    unm_params = []
    for p in params
        if p.managed
            push!(unm_params, param(p.name, ptr_t(p.replacement), false, nothing))
        elseif typeof(p.type) == struct_ref_t
            push!(unm_params, param(p.name, ptr_t(p.type), false, nothing))
        else
            push!(unm_params, p)
        end
    end

    for (ix, p) in enumerate(unm_params)

        cstr *= "\t" * cstring_from_type(p.type) * " " * string(p.name)
        if ix < length(unm_params)
            cstr *= ",\n"
        end
    end
    cstr *= ");"

    return cstr
end

max_stars_repo_name: nec-research/jCall

max_issues_repo_name: nec-research/jCall

max_forks_repo_name: nec-research/jCall


text:
@doc raw"""
```
img = generate_surface(shape::AbstractSyntheticShape, albedo::Real = 1.0, illumination_direction::Vector{T} where T <: Real = [0.0, 0.0, 1.0]; img_size::Int=151, noise::Real = 0, relative_noise::Bool = true)
```
Creates a synthetic image of a given shape with the given illumination properties.
# Output
Returns a grayscale image with dimensions `img_size` in both axes.
# Details
Uses the passed properties to calculate the partial differentials of the shape
reflectance and uses these to generate a synthetic image of the shape to be used
as a consistent test image for SFS testing.
# Arguments
The function arguments are described in more detail below.
## `shape`
A `AbstractSyntheticShape` that specifies the type of shape to be generated. Available
shapes are; SynthSphere, Prism, Ripple, Ripple2, Cake, Cake2, SynthVase, Tent, Dem,
SynthGaussian and SuperGaussian. Each has a parameter `radius`, a `Real` that
specifies the radius of the desired shape or its scale within the image. This can
acts as either an actual radius (such as for `SynthSphere`) or rather as a scale
factor (such as for `Ripple`). If left unspecified a default value is used for each.
## `albedo`
A `Real` that specifies the albedo (amount of light reflected) of the image. If
left unspecified a default value of 1.0 is used.
## `illumination_direction`
A `Vector{T} where T <: Real` that specifies the tilt value to be used by the
algorithm. The `illumination_direction` should be a vector of the form [x,y,z]
where x,y,z are in the range [0,1]. If left unspecified a default value of
[0.0,0.0,1.0] is used.
## `img_size`
A `Int` that specifies the dimentions of the final image should be. Defults to 151.
## `noise`
A `Real` that specifies the standard deviation of noise to be added to the image.
Defults to 0.0.
## `relative_noise`
A `Bool` that specifies if the noise added to the image should be relative to the
gradient of the object displayed in the image at that point.
# Example
Compute the heightmap for a synthetic image generated by `generate_surface`.
```julia
using Images, ShapeFromShading

#generate synthetic image
img = generate_surface(SynthSphere(50), 1.0, [0.5,0.1,0.7])
```
# Reference
1. S. Elhabian, "Hands on Shape from Shading", Computer Vision and Image Processing, 2008.
"""
function generate_surface(shape::AbstractSyntheticShape, albedo::Real = 1.0, illumination_direction::Vector{T} where T <: Real = [0.0, 0.0, 1.0]; img_size::Int=151, noise::Real = 0, relative_noise::Bool = true)
    I = illumination_direction./norm(illumination_direction)
    ρ = albedo
    R = zeros(Float64, img_size, img_size)
    x,y = setup_xy(img_size)
    p,q = synthetic_gradient(shape, img_size=img_size, noise=noise, relative_noise=relative_noise)
    for i in CartesianIndices(R)
        R[i] = (ρ * (-I[1] * p[i] - I[2] * q[i] + I[3])) / sqrt(1 + p[i]^2 + q[i]^2)
    end
    E = max.(0.0, Float64.(R))

    #convert to img and return
    img = Gray.(E)
    return img
end

# Generates image from given gradients
function generate_surface(p, q, albedo::Real = 0.5, illumination_direction::Vector{T} where T <: Real = [0, 0, 1]; img_size::Int=151, noise::Real = 0, relative_noise::Bool = true)
    I = illumination_direction./norm(illumination_direction)
    ρ = albedo
    R = zeros(Float64, img_size, img_size)
    x,y = setup_xy(img_size)
    for i in CartesianIndices(R)
        R[i] = (ρ * (-I[1] * p[i] - I[2] * q[i] + I[3])) / sqrt(1 + p[i]^2 + q[i]^2)
    end
    E = max.(0.0, Float64.(R))

    #convert to img and return
    img = Gray.(E)
    return img
end

# Generates image, gradients and ground truth for a shape
function generate_surface_data(shape::AbstractSyntheticShape, albedo::Real = 0.5, illumination_direction::Vector{T} where T <: Real = [0, 0, 1]; img_size::Int=151, noise::Real = 0, relative_noise::Bool = true)
    img = generate_surface(shape, albedo, illumination_direction, img_size = img_size, noise = noise, relative_noise = relative_noise)
    p,q = synthetic_gradient(shape, img_size = img_size, noise = noise, relative_noise = relative_noise)
    Z = ground_truth(shape, img_size = img_size)
    return img, p, q, Z
end

# Generate three sythetic images of a shape
function generate_photometric(shape::AbstractSyntheticShape, albedo::Real = 0.5, I₁::Vector{T} where T <: Real = [0, 0, 1], I₂::Vector{T} where T <: Real = [0.5, 0, 1], I₃::Vector{T} where T <: Real = [0, 0.5, 1]; radius::Real = 50, img_size::Int=151, noise::Real = 0, relative_noise::Bool = false)
    img1 = generate_surface(shape, albedo, I₁, img_size = img_size, noise = noise, relative_noise = relative_noise)
    img2 = generate_surface(shape, albedo, I₂, img_size = img_size, noise = noise, relative_noise = relative_noise)
    img3 = generate_surface(shape, albedo, I₃, img_size = img_size, noise = noise, relative_noise = relative_noise)
    return img1, img2, img3
end

# Generates synthetic gradients for a shape
function synthetic_gradient(shape::AbstractSyntheticShape; img_size::Int=151, noise::Real = 0, relative_noise::Bool = true)
    od = OnceDifferentiable(x -> shape(x), [0.0,0.0]; autodiff = :forward)
    p = zeros(Float64, img_size, img_size)
    q = zeros(Float64, img_size, img_size)
    x,y = setup_xy(img_size)
    for i in CartesianIndices(x)
        pq = [0.0,0.0]
        od.df(pq,([x[i],y[i]]))
        p[i], q[i] = pq
        if isnan(p[i])
            p[i] = 0.0
        end
        if isnan(q[i])
            q[i] = 0.0
        end
    end
    valP = maximum(abs.(p))
    valQ = maximum(abs.(q))
    for i in CartesianIndices(x)
        if relative_noise == false
            p[i] += rand(Normal(0, noise))
            q[i] += rand(Normal(0, noise))
        else
            p[i] += rand(Normal(0, noise*valP))
            q[i] += rand(Normal(0, noise*valQ))
        end
    end
    return p,q
end

# Generates gorund truth for a shape
function ground_truth(shape::AbstractSyntheticShape;img_size::Int=151)
    Z = zeros(Float64, img_size, img_size)
    x,y = setup_xy(img_size)
    for i in CartesianIndices(Z)
        Z[i] = shape([x[i],y[i]])
    end
    return Z
end

# Cake like shape
function (shape::Cake)(x::AbstractArray)
    r = shape.radius
    z = 75 * sin(π*(x[1]^2 + x[2]^2) / (71 * r))
    if z < 0.0
        z = 0.0
    end
    if x[1]^2 + x[2]^2 > r^2
        z = 0.0
    end
    return z
end

# Cake like shape version 2
function (shape::Cake2)(x::AbstractArray)
    r = shape.radius
    xy = sqrt(x[1]^2 + x[2]^2)
    z = -(xy + 60)*(xy + 20)*(xy - 20)*(xy - 60) / 60000
    if z < 0.0
        z = 0.0
    elseif z > 30.0
        z = 30.0
    end
    return z
end

function g(x,y)
    return 3 * (1 - x)^2 * exp(-x^2 - (y + 1)^2) - 10 * (x / 5 - x^3 - y^5) * exp(-x^2 - y^2) - (1/3) * exp(-(x + 1)^2 - y^2)
end

# DEM
function (shape::Dem)(x::AbstractArray)
    r = shape.radius
    x = x ./ r .* 6.4
    x = x ./ 1.6
    z = g(x[1], x[2]) *r/6.4
    return z
end

# Prism like shape
function (shape::Prism)(x::AbstractArray)
    r = shape.radius
    z = 0.0
    if abs(x[1]) <= r/2 && abs(x[2]) <= r/2
        z = x[1]
    end
    return z
end

# Ripple
function (shape::Ripple)(x::AbstractArray)
    r = shape.radius
    z = 25*cos(sqrt(x[1]^2 + x[2]^2) / r)
    return z
end

# Ripple with decayed edges
function (shape::Ripple2)(x::AbstractArray)
    r = shape.radius
    z = 25 * cos(sqrt(x[1]^2 + x[2]^2)/r) * 62.5 * (1/sqrt(pi * 1250)) * exp(-((x[1]^2 + x[2]^2) / 1250))
    return z
end

# Gaussian to a higher power
function (shape::SuperGaussian)(x::AbstractArray)
    r = shape.radius
    x = x ./ r .* 4
    P = 2
    z = exp(-((x[1]^2)/2+(x[2]^2)/2)^P)*r/4
    return z
end

# Gaussian
function (shape::SynthGaussian)(x::AbstractArray)
    r = shape.radius
    x = x ./ r .* 4
    z = exp(-((x[1]^2)/2+(x[2]^2)/2)) * r/4
    return z
end

# Sphere
function (shape::SynthSphere)(x::AbstractArray)
    r = shape.radius
    if x[1]^2 + x[2]^2 < r^2
        z = sqrt(r^2-(x[1]^2 + x[2]^2))
    else
        z = 0.0
    end
    return z
end

function f(x)
    return -138.24 * x^6 + 92.16 * x^5 + 84.48 * x^4 - 48.64 * x^3 - 17.60 * x^2 + 6.4 * x + 3.2
end

# Vase
function (shape::SynthVase)(x::AbstractArray)
    r = shape.radius
    x = x./r .* 6.4
    if f(x[1]/12.8)^2 >= x[2]^2
        z = sqrt(f(x[1]/12.8)^2 - x[2]^2)*r/6.4
    else
        z = 0.0
    end
    return z
end

# Canadian tent
function (shape::Tent)(x::AbstractArray)
    r = shape.radius
    x = x ./ r .* 6.4
    if abs(x[1]) <= 5.12 && abs(x[2]) <= 5.12
        z = min(-2 * abs(x[1]) + 10.24, -abs(x[2]) + 5.12)*r/6.4
    else
        z = 0.0
    end
    return z
end

# Allows for interactive generation of images
function create_image()
    slider1 = widget(-1:0.01:1)
    slider2 = widget(-1:0.01:1)
    slider3 = widget(0.01:0.01:1)
    slider4 = widget(1:100, value = 5)
    sizeVal = spinbox(1:5000, label="Image size:", value = 151)
    r = slider4[]
    I = [slider1[],slider2[],slider3[]]
    imgSize = sizeVal[]
    img = generate_surface(SynthSphere(), 1, I, radius = 5)
    dis = Observable(img)
    shape = SynthSphere()
    shapeOptions = dropdown(OrderedDict("Sphere"=>SynthSphere(), "Ripple"=>Ripple(), "Ripple2"=>Ripple2(), "Cake"=>Cake(), "Cake2"=>Cake2(), "Vase"=>SynthVase(), "Dem"=>Dem(), "Tent"=>Tent(), "Gaussian"=>SynthGaussian(), "Super Gaussian"=>SuperGaussian()))
    ui = Interact.vbox(slider1,slider2,slider3, dis, Interact.hbox(shapeOptions, sizeVal),slider4)
    w = Window()
    body!(w,ui)

    on(slider1) do val
        I = [slider1[],slider2[],slider3[]]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end

    on(slider2) do val
        I = [slider1[],slider2[],slider3[]]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end

    on(slider3) do val
        I = [slider1[],slider2[],slider3[]]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end

    on(slider4) do val
        r = slider4[]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end

    on(shapeOptions) do val
        shape = shapeOptions[]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end

    on(sizeVal) do val
        imgSize = sizeVal[]
        img = generate_surface(shape, 1.0, I, radius = r, img_size = imgSize)
        dis[] = img
    end
end

ext: jl


text:
staticroot = "/static"  #no trailing slash
site_url = "http://example.com" # for generating absolute urls (eg. sitemap, rss)


text:
# Use baremodule to shave off a few KB from the serialized `.ji` file
baremodule CUDA_compat_jll
using Base
using Base: UUID
using LazyArtifacts
import JLLWrappers

JLLWrappers.@generate_main_file_header("CUDA_compat")
JLLWrappers.@generate_main_file("CUDA_compat", UUID("340b3129-3577-50c6-bc43-5dfc29cb9805"))
end  # module CUDA_compat_jll


text:
function convert(::Type{T}, num::ArbRadixFloat{radix,precision,Tdigit}) where {T<:Number,radix,precision,Tdigit}
    y = zero(T)
    monomial = radix^num.exponent
    for digit in num.significand
        y += digit * monomial
        monomial /= radix
    end
    y
end

convert(::Type{ArbRadixFloat{radix,precision}}, num :: Number, verbose::Bool=false) where {radix,precision} =
    convert(ArbRadixFloat{radix,precision,Int}, num)

convert(::Type{ArbRadixFloat{radix,precision,Tdigit}}, num :: ArbRadixFloat{radix,precision,Tdigit},
        verbose::Bool=false) where {radix,precision,Tdigit} = num

function convert(::Type{ArbRadixFloat{radix,precision,Tdigit}},
                      num :: Number) where {radix,precision,Tdigit}
    if abs(radix) >= 1
        exponent = find_exponent(num, radix)
        significand = factorize_leastsquares(num, radix, precision, Tdigit, exponent)
        ArbRadixFloat{radix,precision,Tdigit}(significand, exponent)
    else #Compute using inverse radix and then convert back
        iradix = inv(radix)
        exponent = find_exponent(num, iradix)
        significand = factorize_leastsquares(num, iradix, precision, Tdigit, exponent)
        ArbRadixFloat{radix,precision,Tdigit}(reverse!(significand), precision-1-exponent)
    end
end

##############################################################################
# Unexported private methods for conversion
##############################################################################

function find_exponent(num, radix)
    log_abs_radix = log(abs(radix))
    if log_abs_radix == 0
        #Do not attempt to normalize significand if radix has magnitude 1
        exponent = 0
    else #log_abs_radix != 0
        exponent = floor(Int, log(abs(num))/log_abs_radix)
    end
    exponent
end

function factorize_leastsquares(num, radix, precision, Tdigit, exponent;
                   verbose::Bool=false)

    significand = zeros(Tdigit, precision)

    # Super expensive greedy exhaustion following a full least squares solve
    # This is completely overkill because
    # D is low rank and symmetric (rank 1 for real radix, 2 for complex)
    # we only need the first entry of signif
    # and also "wrong" because the least squares solution is "undesirable"
    # because if you look at unit radixes like im, you will see the solution
    # norm spread out over many elements but really you just want the sparsest one
    # for unit radixes there is usually no unique numeration
    # and this method will eventually converge on one particular numeration
    # with the mass focused toward the last few 'digits'
    # but in some sense it's not the intuitive one; want the mass focused toward
    # the front
    if verbose
        println("Starting norm: ", abs(num))
        println("Starting exponent: ", exponent)
        println("Power\t| Digit\t| Solution \t| Residual norm\t| Change")
        println("-----\t| -----\t| -------- \t| -------------\t| ------")
    end

    #Element type of working intermediate problem
    eltype = typeof(float(real(radix)))

    for digit_id=1:precision
        nk = precision+1-digit_id
        D = zeros(eltype, nk, nk)
        c = zeros(eltype, nk)
        r = float(radix)
        ir = inv(r)
        irconj = conj(ir)
        for i=1:nk, j=1:nk
            z = irconj^(i+digit_id-1)*ir^(j+digit_id-1)
            D[i,j] = real(z+conj(z))
        end
        for i=1:nk
            z = irconj^(i+digit_id-1)*num*ir^(exponent+1)
            c[i] = real(z+conj(z))
        end
        signif = pinv(D)*c
        significand[digit_id] = digit = round(Tdigit, signif[1])
        powr = exponent+1-digit_id
        if powr >= 0
            Δ = digit * r^powr
        else #Some number types don't have x^-n defined, manually rewrite as inv(x)^n
            Δ = digit * ir^-powr
        end
        num -= Δ
        if verbose
            println(powr, "\t|", digit, "\t|", signif, "\t|", abs(num), "\t|", abs(Δ))
        end
        if num == 0 #Exactly representable
            if verbose
                println("Exact representation found")
            end
            break
        end
    end
    significand
end

function factorize_greedy(num, radix, precision, Tdigit, exponent;
                   verbose::Bool=false)

    significand = zeros(Tdigit, precision)
    monomial = radix ^ exponent
    if verbose
        println("Starting norm: ", abs(num))
        println("Starting exponent: ", exponent)
        println("Power\t| Digit\t| Quotient\t| Residual norm\t| Change")
        println("-----\t| -----\t| --------\t| -------------\t| ------")
    end
    for digit_id in 1:precision
        quotient = real(conj(num)*monomial)/real(conj(monomial)*monomial)
        digit = round(Tdigit, quotient)
        num -= digit*monomial
        monomial /= radix
        significand[digit_id] = digit
        if verbose
            println(exponent+1-digit_id, "\t|", digit, "\t|", quotient, "\t|", abs(num), "\t|", abs(digit*monomial))
        end
        if num == 0 #Exactly representable
            if verbose
                println("Exact representation found")
            end
            break
        end
    end
    significand
end


text:
module GridVisualize

using Printf
using LinearAlgebra

using DocStringExtensions
using OrderedCollections
using ElasticArrays
using StaticArrays
using Colors
using ColorSchemes
using GeometryBasics
using PkgVersion

using ExtendableGrids
using Requires

function __init__()
    # Set default backend depending on installed packages
    global default_backend=nothing
    @require Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" default_backend=Plots
    @require PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee" default_backend=PyPlot
    @require GLMakie = "e9467ef8-e4e7-5192-8a1a-b1aee30e663a" default_backend=GLMakie
    @require PlutoVista = "646e1f28-b900-46d7-9d87-d554eb38a413" default_backend=PlutoVista
end

include("dispatch.jl")
include("common.jl")
include("pyplot.jl")
include("makie.jl")
include("vtkview.jl")
include("meshcat.jl")
include("plots.jl")
include("plutovista.jl")

export scalarplot,scalarplot!
export gridplot,gridplot!
export vectorplot,vectorplot!
export streamplot,streamplot!
export save,reveal,backend!
export isplots,isvtkview,ispyplot,ismakie,isplutovista
export GridVisualizer, SubVisualizer
export plottertype, available_kwargs
export default_plotter!, default_plotter
export PyPlotType,MakieType,PlotsType,VTKViewType, PlutoVistaType, MeshCatType

end


text:
# Test of binary heaps

using DataStructures
using Base.Test

# auxiliary functions

function heap_values{VT,Comp}(h::MutableBinaryHeap{VT,Comp})
    n = length(h)
    nodes = h.nodes
    @assert length(nodes) == n
    vs = Array(VT, n)
    for i = 1 : n
        vs[i] = nodes[i].value
    end
    vs
end

function list_values{VT,Comp}(h::MutableBinaryHeap{VT,Comp})
    n = length(h)
    nodes = h.nodes
    nodemap = h.node_map
    vs = Array(VT, 0)
    for i = 1 : length(nodemap)
        id = nodemap[i]
        if id > 0
            push!(vs, nodes[id].value)
        end
    end
    vs
end

function verify_heap{VT,Comp}(h::MutableBinaryHeap{VT,Comp})
    comp = h.comparer
    nodes = h.nodes
    n = length(h)
    m = n / 2
    for i = 1 : m
        v = nodes[i].value
        lc = i * 2
        if lc <= n
            if compare(comp, nodes[lc].value, v)
                return false
            end
        end
        rc = lc + 1
        if rc <= n
            if compare(comp, nodes[rc].value, v)
                return false
            end
        end
    end
    return true
end

# test make heap

vs = [4, 1, 3, 2, 16, 9, 10, 14, 8, 7]
h = mutable_binary_minheap(vs)

@test length(h) == 10
@test !isempty(h)
@test top(h) == 1
@test isequal(list_values(h), vs)
@test isequal(heap_values(h), [1, 2, 3, 4, 7, 9, 10, 14, 8, 16])

h = mutable_binary_maxheap(vs)

@test length(h) == 10
@test !isempty(h)
@test top(h) == 16
@test isequal(list_values(h), vs)
@test isequal(heap_values(h), [16, 14, 10, 8, 7, 3, 9, 1, 4, 2])

# test push!

hmin = mutable_binary_minheap(Int)
@test length(hmin) == 0
@test isempty(hmin)

ss = {
    [4],
    [1, 4],
    [1, 4, 3],
    [1, 2, 3, 4],
    [1, 2, 3, 4, 16],
    [1, 2, 3, 4, 16, 9],
    [1, 2, 3, 4, 16, 9, 10],
    [1, 2, 3, 4, 16, 9, 10, 14],
    [1, 2, 3, 4, 16, 9, 10, 14, 8],
    [1, 2, 3, 4, 7, 9, 10, 14, 8, 16]}

for i = 1 : length(vs)
    ia = push!(hmin, vs[i])
    @test ia == i
    @test length(hmin) == i
    @test !isempty(hmin)
    @test isequal(list_values(hmin), vs[1:i])
    @test isequal(heap_values(hmin), ss[i])
end

hmax = mutable_binary_maxheap(Int)
@test length(hmax) == 0
@test isempty(hmax)

ss = {
    [4],
    [4, 1],
    [4, 1, 3],
    [4, 2, 3, 1],
    [16, 4, 3, 1, 2],
    [16, 4, 9, 1, 2, 3],
    [16, 4, 10, 1, 2, 3, 9],
    [16, 14, 10, 4, 2, 3, 9, 1],
    [16, 14, 10, 8, 2, 3, 9, 1, 4],
    [16, 14, 10, 8, 7, 3, 9, 1, 4, 2]}

for i = 1 : length(vs)
    ia = push!(hmax, vs[i])
    @test ia == i
    @test length(hmax) == i
    @test !isempty(hmax)
    @test isequal(list_values(hmax), vs[1:i])
    @test isequal(heap_values(hmax), ss[i])
end

# test pop!

@test isequal(extract_all!(hmin), [1, 2, 3, 4, 7, 8, 9, 10, 14, 16])
@test isempty(hmin)

@test isequal(extract_all!(hmax), [16, 14, 10, 9, 8, 7, 4, 3, 2, 1])
@test isempty(hmax)

# test hybrid push! and pop!

h = mutable_binary_minheap(Int)

push!(h, 5)
push!(h, 10)
@test isequal(heap_values(h), [5, 10])
@test isequal(list_values(h), [5, 10])

@test pop!(h) == 5
@test isequal(heap_values(h), [10])
@test isequal(list_values(h), [10])

push!(h, 7)
push!(h, 2)
@test isequal(heap_values(h), [2, 10, 7])
@test isequal(list_values(h), [10, 7, 2])

@test pop!(h) == 2
@test isequal(heap_values(h), [7, 10])
@test isequal(list_values(h), [10, 7])

# test update!

xs = rand(100)
hmin = mutable_binary_minheap(xs)
hmax = mutable_binary_maxheap(xs)

@test length(hmin) == 100
@test length(hmax) == 100
@test verify_heap(hmin)
@test verify_heap(hmax)

for t = 1 : 1000
    i = rand(1:100)
    v = rand()
    xs[i] = v
    update!(hmin, i, v)
    update!(hmax, i, v)
    @test length(hmin) == 100
    @test length(hmax) == 100
    @test verify_heap(hmin)
    @test verify_heap(hmax)
    @test isequal(list_values(hmin), xs)
    @test isequal(list_values(hmax), xs)
end


text:
"""
    lambdipole(U, R, g::TwoDimGrid; center=(mean(g.x), mean(g.y))

Return the 2D vorticity field of the Lamb dipole with strength `U` and radius `R`, centered on
`center=(xc, yc)` and on the grid `g`. The default value of `center` is the middle of the grid.
"""
function lambdipole(U, R, g::AbstractGrid{T, A}; center=(mean(g.x), mean(g.y))) where {T, A}
  firstzero = 3.8317059702075123156
  k = firstzero/R # dipole wavenumber for radius R in terms of first zero of besselj
  q0 = -2U*k/besselj(0, k*R) # dipole amplitude for strength U and radius R

  xc, yc = center
  r = @. sqrt( (g.x-xc)^2 + (g.y-yc)^2 )
  besselj1 = A([besselj(1, k*r[i, j]) for i=1:g.nx, j=1:g.ny])
  q = @. q0 * besselj1 * (g.y-yc)/r
  @. q[r >= R] = 0
  q
end

"""
    peakedisotropicspectrum(g, kpeak, E0; mask=mask, allones=false)
Generate a real and random two-dimensional vorticity field q(x, y) with
a Fourier spectrum peaked around a central non-dimensional wavenumber kpeak and
normalized so that its total energy is E0.
"""
function peakedisotropicspectrum(g::AbstractGrid{T, A}, kpeak::Real, E0::Real; mask=ones(size(g.Krsq)), allones=false) where {T, A}
  if g.Lx !== g.Ly
    error("the domain is not square")
  else
    k0 = kpeak*2π/g.Lx
    modk = sqrt.(g.Krsq)
    psik = A(zeros(T, (g.nk, g.nl)))
    psik = @. (modk^2 * (1 + (modk/k0)^4))^(-0.5)
    psik[1, 1] = 0.0
    phases = randn(Complex{T}, size(g.Krsq))
    phases_real, phases_imag = real.(phases), imag.(phases)
    phases = A(phases_real) + im*A(phases_imag)
    psih = @. phases*psik
    if allones; psih = psik; end
    psih = psih.*A(mask)
    Ein = real(sum(g.Krsq.*abs2.(psih)/(g.nx*g.ny)^2))
    psih = psih*sqrt(E0/Ein)
    q = A(-irfft(g.Krsq.*psih, g.nx))
  end
end


text:
test = """939
7,13,x,x,59,x,31,19
"""

function parse_input(input::AbstractString)::Tuple
    (departure, buses) = split(input, r"\s+"; keepempty = false)
    buses = tryparse.(Int, split(buses, ','))
    return (parse(Int, departure), buses)
end

function q1(departure::Integer, buses::AbstractVector)::Int
    buses = filter(!isnothing, buses)
    bus = buses[argmin([bus - departure % bus for bus in buses])]
    return (bus - departure % bus) * bus
end

function q2(buses::AbstractVector)::Int
    return mod_lcm([(bus, rem(1 - i, bus, RoundDown))
                    for (i, bus) in enumerate(buses)
                    if bus !== nothing])
end

function mod_lcm(vs::AbstractVector)::Int
    isempty(vs) && return 0
    ((v1, m1), rest...) = vs
    r = map(rest) do (v, m)
        vp = v ÷ gcd(v, v1)
        k = findfirst(==(m), [(k * v1 + m1) % v for k in 0:vp - 1]) - 1
        return (vp, k)
    end |> mod_lcm
    return r * v1 + m1
end

let
    (departure, buses) = parse_input(test)
    @assert q1(departure, buses) == 295
    @assert q2(buses) == 1068781
end

(departure, buses) = parse_input(read("day13.in", String))
println("Q1: ", q1(departure, buses))
println("Q1: ", q2(buses))

max_stars_repo_name: dhanak/competitive-coding

max_issues_repo_name: dhanak/competitive-coding

max_forks_repo_name: dhanak/competitive-coding


text:
@testset "get_bdedge_integration_pts" begin
    mesh = Mesh(1, 1, 1.0)
    bd = bcedge(mesh)
    pts = get_bdedge_integration_pts(bd, mesh)
    @test maximum(abs.(pts - [ 0.788675  0.0
            0.211325  0.0
            1.0       0.788675
            1.0       0.211325
            0.0       0.788675
            0.0       0.211325
            0.788675  1.0
            0.211325  1.0])) < 1e-5
end

@testset "traction" begin

    for d in [1, 2]
        n = 20
        mesh = Mesh(n, n, 1/n; degree=d)
        bd = bcedge(mesh)

        # lower
        bdedge = []
        for i = 1:size(bd, 1)
            a, b = bd[i,:]
            if mesh.nodes[a,2]<1e-5 && mesh.nodes[b, 2]<1e-5
                push!(bdedge, [a b])
            end
        end
        bdedge = vcat(bdedge...)
        t = eval_f_on_boundary_edge((x, y)->y^2 + x^2 + 1.0, bdedge, mesh)
        rhs = compute_fem_traction_term1(t, bdedge, mesh)
        @test abs(sum(rhs)-4/3)<1e-5

        # upper
        bdedge = []
        for i = 1:size(bd, 1)
            a, b = bd[i,:]
            if mesh.nodes[a,2]>1-1e-5 && mesh.nodes[b, 2]>1-1e-5
                push!(bdedge, [a b])
            end
        end
        bdedge = vcat(bdedge...)
        t = eval_f_on_boundary_edge((x, y)-> y^2 + x^2 + 1.0, bdedge, mesh)
        rhs = compute_fem_traction_term1(t, bdedge, mesh; order=6)
        @test abs(sum(rhs)-7/3)<1e-5

        # left
        bdedge = []
        for i = 1:size(bd, 1)
            a, b = bd[i,:]
            if mesh.nodes[a,1]<1e-5 && mesh.nodes[b, 1]<1e-5
                push!(bdedge, [a b])
            end
        end
        bdedge = vcat(bdedge...)
        t = eval_f_on_boundary_edge((x, y)->y^2 + x^2 + 1.0, bdedge, mesh)
        rhs = compute_fem_traction_term1(t, bdedge, mesh)
        @test abs(sum(rhs)-4/3)<1e-5

        # right
        bdedge = []
        for i = 1:size(bd, 1)
            a, b = bd[i,:]
            if mesh.nodes[a,1]>1-1e-5 && mesh.nodes[b, 1]>1-1e-5
                push!(bdedge, [a b])
            end
        end
        bdedge = vcat(bdedge...)
        t = eval_f_on_boundary_edge((x, y)->y^2 + x^2 + 1.0, bdedge, mesh)
        rhs = compute_fem_traction_term1(t, bdedge, mesh)
        @test abs(sum(rhs)-7/3)<1e-5
    end
end

@testset "impose_Dirichlet_boundary_conditions" begin
    A = rand(5, 5)
    bd = [5;2;3]
    bdval = Float64[1;2;3]
    rhs = rand(5)
    A1, rhs1 = impose_Dirichlet_boundary_conditions(A, rhs, bd, bdval)
    A2, rhs2 = impose_Dirichlet_boundary_conditions(constant(sparse(A)), rhs, bd, bdval)

    @test maximum(abs.(Array(run(sess, A2)) - A1)) < 1e-5
    @test maximum(abs.(run(sess, rhs2) - rhs1)) < 1e-5
end


text:
export TimeDomainData, ComplexSpectrum

struct TimeDomainData <: TimeResolvedData
    td::Vector{Float64} # time delay
    ΔE::Vector{Float64} # electric field difference
    Δt::Float64

    function TimeDomainData(fname::String; args...)
        dat = readdlm(fname; args...)
        td  = collect(range(dat[1,1], dat[end,1], length=clp2(size(dat, 1))))
        ΔE  = similar(td)
        polyinterp!(ΔE, td, view(dat, :, 2), view(dat, :, 1))
        return new(td, ΔE, td[2] - td[1])
    end
end

struct ComplexSpectrum <: TimeResolvedData
    @def prop Vector{Float64} fs Re Im Rs θs; sz::Int

    function ComplexSpectrum(ref::NamedTuple{(:frq, :cxy), Tuple{Vector{Float64}, Matrix{Float64}}},
                             sam::NamedTuple{(:frq, :cxy), Tuple{Vector{Float64}, Matrix{Float64}}};
                             fmin::Real=0.1, fmax::Real=3.0)
        frq = ref.frq
        length(frq) ≠ length(sam.frq) && error("length(ref.frq) ≠ length(sam.frq)")

        lx = length(frq) >> 1
        @inbounds while frq[lx] < fmin
            lx += 1
        end
        rx = lx
        @inbounds while !(frq[rx] > fmax)
            rx += 1
        end
        sz = rx - lx

        @def vars Vector{Float64}(undef, sz) fs Re Im Rs θs

        @simd for i in eachindex(fs)
            @inbounds fs[i] = frq[lx+i]
        end

        rxy = ref.cxy
        sxy = sam.cxy

        @inbounds for i in eachindex(fs)
            Re[i], Im[i] = cdiv(sxy[lx+i,1], sxy[lx+i,2], rxy[lx+i,1], rxy[lx+i,2])
        end

        #### Phases unwrapping
        halfπ = 0.5 * π
        @inbounds θs[1] = prevθ = atan(Im[1], Re[1])
        @inbounds Rs[1] = 0.0
        @inbounds for i in 2:sz
            θs[i] = thisθ = atan(Im[i], Re[i])
            diff = thisθ - prevθ
            θmod = rem2pi(diff + halfπ, RoundDown) - halfπ
            diff > 0.0 && θmod == -halfπ && (θmod = halfπ)
            Rs[i] = abs(diff) < halfπ ? Rs[i-1] : θmod - diff + Rs[i-1]
            prevθ = thisθ
        end

        @simd for i in eachindex(fs)
            @inbounds θs[i] += Rs[i]
        end

        @simd for i in eachindex(fs)
            @inbounds Rs[i] = apy2(Re[i], Im[i])
        end

        return new(fs, Re, Im, Rs, θs, sz)
    end
end


text:
module SocialSystems

using StatsBase: weights, sample
using Distributions: MvNormal, rand, Normal, cdf, pdf

import Base: +, -, *, length, size, getindex, setindex!, start, done, next, show, dot

# Utils Export
export randSphere, gammasoc, phi, G

# Types Export
export MoralVector, Society, BasicSociety, DistrustSociety
export agents, interactions, insertagent!, hamiltonian, magnetization, believeness, quadrupole, consensus
export epssoc, rhosoc, gammasoc, cogcost

# Dynamics Export
export metropolisStep!, metropolis!,
        discreteStep!, discreteEvol!, societyHistory!,
        computeDeltas

include("constants.jl")

include("utils.jl")

include("types/MoralVector.jl")
include("types/Society.jl")
include("types/BasicSociety.jl")
include("types/DistrustSociety.jl")

include("dynamics.jl")

end #module


text:
mutable struct Priors
  delta0
  delta1
  sig2
  W
  eta0
  eta1
  v
  H
  alpha
  eps
end

# NOTE: Uses Gamma(shape, scale) parameterization.

function Priors(;K::Int, L::Dict{Bool, Int}, use_stickbreak::Bool=false)
  if use_stickbreak
    v_prior = a -> Beta(a, oftype(a, 1))
  else
    v_prior = a -> Beta(a / K, oftype(a, 1))
  end

  return Priors(Gamma(1, 1),# delta0
                Gamma(1, 1), # delta1
                LogNormal(-1, .1), # sig2
                # Gamma(2, .1), # sig2
                Dirichlet(ones(K) / K), # W
                Dirichlet(ones(L[0]) / L[0]), # eta0
                Dirichlet(ones(L[1]) / L[1]), # eta1
                v_prior, # v
                Uniform(0, 1), # H
                Gamma(0.1, 10.0), # alpha
                Beta(1, 99)) # eps
end

max_stars_repo_name: luiarthur/CytofResearch

max_issues_repo_name: luiarthur/CytofResearch

max_forks_repo_name: luiarthur/CytofResearch


text:
abstract type ParentLinks end

"""
  Indicates that this tree stores parent links explicitly. The implementation
  is responsible for defining the parentind function to expose this
  information.
"""
struct StoredParents <: ParentLinks; end
struct ImplicitParents <: ParentLinks; end

parentlinks(::Type) = ImplicitParents()
parentlinks(tree) = parentlinks(typeof(tree))

abstract type SiblingLinks end

"""
  Indicates that this tree stores sibling links explicitly, or can compute them
  quickly (e.g. because the tree has a (small) fixed branching ratio, so the
  current index of a node can be determined by quick linear search). The
  implementation is responsible for defining the relative_state function
  to expose this information.
"""
struct StoredSiblings <: SiblingLinks; end
struct ImplicitSiblings <: SiblingLinks; end

siblinglinks(::Type) = ImplicitSiblings()
siblinglinks(tree) = siblinglinks(typeof(tree))

abstract type TreeKind end
struct RegularTree <: TreeKind; end
struct IndexedTree <: TreeKind; end

treekind(tree::Type) = RegularTree()
treekind(tree) = treekind(typeof(tree))
children(tree, node, ::RegularTree) = children(node)
children(tree, node, ::IndexedTree) = (tree[y] for y in childindices(tree, node))
children(tree, node) = children(tree, node, treekind(tree))

function rootstate()
end

childindices(tree, node) =
  tree == node ? childindices(tree, rootstate(tree)) :
  error("Must implement childindices(tree, node)")
function childindices()
end

function parentind
end

"""
    nodetype(tree)

A trait function, defined on the tree object, specifying the types of the nodes.
The default is `Any`. When applicable, define this trait to make iteration inferrable.

# Example

```jldoctest
struct IntTree
    num::Int
    children::Vector{IntTree}
end
AbstractTrees.children(itree::IntTree) = itree.children
AbstractTrees.nodetype(::IntTree) = IntTree
```

This suffices to make iteration over, e.g., `Leaves(itree::IntTree)` inferrable.
"""
nodetype(tree) = Any
idxtype(tree) = Int


text:
export AbstractSystem, AbstractHomotopy, jacobian

"""
    AbstractSystem

An abstract type representing a polynomial system ``F(x, p)`` where ``x`` are variables
and ``p`` are possible parameters.
"""
abstract type AbstractSystem end

# Base.size(F::AbstractSystem)
Base.size(F::AbstractSystem, i::Integer) = size(F)[i]

"""
    parameters(F::AbstractSystem)

Returns the [`Variable`](@ref)s used in the system.
"""
parameters(F::AbstractSystem) = Variable[]

"""
    variable_groups(F::AbstractSystem)

Returns the variable groups of the system.
"""
variable_groups(::AbstractSystem) = nothing

"""
    nvariables(F::AbstractSystem)

Returns the number of variables of a given system `F`.
"""
nvariables(F::AbstractSystem) = size(F, 2)

"""
    nparameters(F::AbstractSystem)

Returns the number of parameters of a given system `F`.
"""
nparameters(F::AbstractSystem) = length(parameters(F))

"""
    System(F::AbstractSystem)

Construct a (symbolic) [`System`](@ref) from `F`.
"""
function System(F::AbstractSystem)
    x, p = variables(F), parameters(F)
    System(F(x, p), x, p, variable_groups(F))
end

(F::AbstractSystem)(x, p = nothing) = evaluate(F, x, p)

"""
    evaluate(F::AbstractSystem, x, p = nothing)

Evaluate the given system.
"""
function evaluate(F::AbstractSystem, x, p = nothing)
    u = Vector{Any}(undef, size(F, 1))
    evaluate!(u, F, x, p)
    ModelKit.to_smallest_eltype(u)
end

"""
    jacobian(F::AbstractSystem, x, p = nothing)

Compute the Jacobian of the given system.
"""
function jacobian(F::AbstractSystem, x, p = nothing)
    u = Vector{Any}(undef, size(F, 1))
    U = Matrix{Any}(undef, size(F))
    evaluate_and_jacobian!(u, U, F, x, p)
    ModelKit.to_smallest_eltype(U)
end

##############
## Homotopy ##
##############

"""
    AbstractHomotopy

An abstract type representing a homotopy ``H(x, t)`` where ``x`` are variables
and ``t`` is space / time.
"""
abstract type AbstractHomotopy end

Base.size(H::AbstractHomotopy, i::Integer) = i == 1 ? first(size(H)) : last(size(H))

(H::AbstractHomotopy)(x, t, p = nothing) = evaluate(H, x, t, p)
function evaluate(H::AbstractHomotopy, x, t, p = nothing)
    U = Vector{Any}(undef, first(size(H)))
    to_smallest_eltype(evaluate!(U, H, x, t, p))
end

function jacobian(H::AbstractHomotopy, x, t, p = nothing)
    n, m = size(H)
    u = Vector{Any}(undef, size(H, 1))
    U = Matrix{Any}(undef, size(H))
    evaluate_and_jacobian!(u, U, H, x, t, p)
    to_smallest_eltype(U)
end


text:
module PyVirtualenvExample

using PyVirtualenv

function __init__()
    PyVirtualenv.activate_pipenv(@__MODULE__)
    @eval using PyCall
end

function demo()
    sys_executable = pyimport("sys")[:executable]
    @info "Using: sys.executable = $sys_executable"
end

end # module

max_stars_repo_name: tkf/PyVirtualenvExample

max_issues_repo_name: tkf/PyVirtualenvExample

max_forks_repo_name: tkf/PyVirtualenvExample