text:
# syntax: proto2
using ProtoBuf
import ProtoBuf.meta

mutable struct HeaderBBox <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function HeaderBBox(; kwargs...)
        obj = new(meta(HeaderBBox), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct HeaderBBox
const __meta_HeaderBBox = Ref{ProtoMeta}()
function meta(::Type{HeaderBBox})
    ProtoBuf.metalock() do
        if !isassigned(__meta_HeaderBBox)
            __meta_HeaderBBox[] = target = ProtoMeta(HeaderBBox)
            req = Symbol[:left,:right,:top,:bottom]
            wtype = Dict(:left => :sint64, :right => :sint64, :top => :sint64, :bottom => :sint64)
            allflds = Pair{Symbol,Union{Type,String}}[:left => Int64, :right => Int64, :top => Int64, :bottom => Int64]
            meta(target, HeaderBBox, allflds, req, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_HeaderBBox[]
    end
end
function Base.getproperty(obj::HeaderBBox, name::Symbol)
    if name === :left
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :right
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :top
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :bottom
        return (obj.__protobuf_jl_internal_values[name])::Int64
    else
        getfield(obj, name)
    end
end

mutable struct HeaderBlock <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function HeaderBlock(; kwargs...)
        obj = new(meta(HeaderBlock), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct HeaderBlock
const __meta_HeaderBlock = Ref{ProtoMeta}()
function meta(::Type{HeaderBlock})
    ProtoBuf.metalock() do
        if !isassigned(__meta_HeaderBlock)
            __meta_HeaderBlock[] = target = ProtoMeta(HeaderBlock)
            fnum = Int[1,4,5,16,17,32,33,34]
            allflds = Pair{Symbol,Union{Type,String}}[:bbox => HeaderBBox, :required_features => Base.Vector{AbstractString}, :optional_features => Base.Vector{AbstractString}, :writingprogram => AbstractString, :source => AbstractString, :osmosis_replication_timestamp => Int64, :osmosis_replication_sequence_number => Int64, :osmosis_replication_base_url => AbstractString]
            meta(target, HeaderBlock, allflds, ProtoBuf.DEF_REQ, fnum, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_HeaderBlock[]
    end
end
function Base.getproperty(obj::HeaderBlock, name::Symbol)
    if name === :bbox
        return (obj.__protobuf_jl_internal_values[name])::HeaderBBox
    elseif name === :required_features
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{AbstractString}
    elseif name === :optional_features
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{AbstractString}
    elseif name === :writingprogram
        return (obj.__protobuf_jl_internal_values[name])::AbstractString
    elseif name === :source
        return (obj.__protobuf_jl_internal_values[name])::AbstractString
    elseif name === :osmosis_replication_timestamp
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :osmosis_replication_sequence_number
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :osmosis_replication_base_url
        return (obj.__protobuf_jl_internal_values[name])::AbstractString
    else
        getfield(obj, name)
    end
end

mutable struct StringTable <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function StringTable(; kwargs...)
        obj = new(meta(StringTable), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct StringTable
const __meta_StringTable = Ref{ProtoMeta}()
function meta(::Type{StringTable})
    ProtoBuf.metalock() do
        if !isassigned(__meta_StringTable)
            __meta_StringTable[] = target = ProtoMeta(StringTable)
            allflds = Pair{Symbol,Union{Type,String}}[:s => Base.Vector{Vector{UInt8}}]
            meta(target, StringTable, allflds, ProtoBuf.DEF_REQ, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_StringTable[]
    end
end
function Base.getproperty(obj::StringTable, name::Symbol)
    if name === :s
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Vector{UInt8}}
    else
        getfield(obj, name)
    end
end

mutable struct Info <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function Info(; kwargs...)
        obj = new(meta(Info), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct Info
const __meta_Info = Ref{ProtoMeta}()
function meta(::Type{Info})
    ProtoBuf.metalock() do
        if !isassigned(__meta_Info)
            __meta_Info[] = target = ProtoMeta(Info)
            val = Dict{Symbol,Any}(:version => -1)
            allflds = Pair{Symbol,Union{Type,String}}[:version => Int32, :timestamp => Int64, :changeset => Int64, :uid => Int32, :user_sid => UInt32, :visible => Bool]
            meta(target, Info, allflds, ProtoBuf.DEF_REQ, ProtoBuf.DEF_FNUM, val, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_Info[]
    end
end
function Base.getproperty(obj::Info, name::Symbol)
    if name === :version
        return (obj.__protobuf_jl_internal_values[name])::Int32
    elseif name === :timestamp
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :changeset
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :uid
        return (obj.__protobuf_jl_internal_values[name])::Int32
    elseif name === :user_sid
        return (obj.__protobuf_jl_internal_values[name])::UInt32
    elseif name === :visible
        return (obj.__protobuf_jl_internal_values[name])::Bool
    else
        getfield(obj, name)
    end
end

mutable struct DenseInfo <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function DenseInfo(; kwargs...)
        obj = new(meta(DenseInfo), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct DenseInfo
const __meta_DenseInfo = Ref{ProtoMeta}()
function meta(::Type{DenseInfo})
    ProtoBuf.metalock() do
        if !isassigned(__meta_DenseInfo)
            __meta_DenseInfo[] = target = ProtoMeta(DenseInfo)
            pack = Symbol[:version,:timestamp,:changeset,:uid,:user_sid,:visible]
            wtype = Dict(:timestamp => :sint64, :changeset => :sint64, :uid => :sint32, :user_sid => :sint32)
            allflds = Pair{Symbol,Union{Type,String}}[:version => Base.Vector{Int32}, :timestamp => Base.Vector{Int64}, :changeset => Base.Vector{Int64}, :uid => Base.Vector{Int32}, :user_sid => Base.Vector{Int32}, :visible => Base.Vector{Bool}]
            meta(target, DenseInfo, allflds, ProtoBuf.DEF_REQ, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, pack, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_DenseInfo[]
    end
end
function Base.getproperty(obj::DenseInfo, name::Symbol)
    if name === :version
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    elseif name === :timestamp
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :changeset
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :uid
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    elseif name === :user_sid
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    elseif name === :visible
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Bool}
    else
        getfield(obj, name)
    end
end

mutable struct ChangeSet <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function ChangeSet(; kwargs...)
        obj = new(meta(ChangeSet), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct ChangeSet
const __meta_ChangeSet = Ref{ProtoMeta}()
function meta(::Type{ChangeSet})
    ProtoBuf.metalock() do
        if !isassigned(__meta_ChangeSet)
            __meta_ChangeSet[] = target = ProtoMeta(ChangeSet)
            req = Symbol[:id]
            allflds = Pair{Symbol,Union{Type,String}}[:id => Int64]
            meta(target, ChangeSet, allflds, req, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_ChangeSet[]
    end
end
function Base.getproperty(obj::ChangeSet, name::Symbol)
    if name === :id
        return (obj.__protobuf_jl_internal_values[name])::Int64
    else
        getfield(obj, name)
    end
end

mutable struct Node <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function Node(; kwargs...)
        obj = new(meta(Node), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct Node
const __meta_Node = Ref{ProtoMeta}()
function meta(::Type{Node})
    ProtoBuf.metalock() do
        if !isassigned(__meta_Node)
            __meta_Node[] = target = ProtoMeta(Node)
            req = Symbol[:id,:lat,:lon]
            fnum = Int[1,2,3,4,8,9]
            pack = Symbol[:keys,:vals]
            wtype = Dict(:id => :sint64, :lat => :sint64, :lon => :sint64)
            allflds = Pair{Symbol,Union{Type,String}}[:id => Int64, :keys => Base.Vector{UInt32}, :vals => Base.Vector{UInt32}, :info => Info, :lat => Int64, :lon => Int64]
            meta(target, Node, allflds, req, fnum, ProtoBuf.DEF_VAL, pack, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_Node[]
    end
end
function Base.getproperty(obj::Node, name::Symbol)
    if name === :id
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :keys
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :vals
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :info
        return (obj.__protobuf_jl_internal_values[name])::Info
    elseif name === :lat
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :lon
        return (obj.__protobuf_jl_internal_values[name])::Int64
    else
        getfield(obj, name)
    end
end

mutable struct DenseNodes <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function DenseNodes(; kwargs...)
        obj = new(meta(DenseNodes), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct DenseNodes
const __meta_DenseNodes = Ref{ProtoMeta}()
function meta(::Type{DenseNodes})
    ProtoBuf.metalock() do
        if !isassigned(__meta_DenseNodes)
            __meta_DenseNodes[] = target = ProtoMeta(DenseNodes)
            fnum = Int[1,5,8,9,10]
            pack = Symbol[:id,:lat,:lon,:keys_vals]
            wtype = Dict(:id => :sint64, :lat => :sint64, :lon => :sint64)
            allflds = Pair{Symbol,Union{Type,String}}[:id => Base.Vector{Int64}, :denseinfo => DenseInfo, :lat => Base.Vector{Int64}, :lon => Base.Vector{Int64}, :keys_vals => Base.Vector{Int32}]
            meta(target, DenseNodes, allflds, ProtoBuf.DEF_REQ, fnum, ProtoBuf.DEF_VAL, pack, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_DenseNodes[]
    end
end
function Base.getproperty(obj::DenseNodes, name::Symbol)
    if name === :id
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :denseinfo
        return (obj.__protobuf_jl_internal_values[name])::DenseInfo
    elseif name === :lat
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :lon
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :keys_vals
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    else
        getfield(obj, name)
    end
end

mutable struct Way <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function Way(; kwargs...)
        obj = new(meta(Way), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct Way
const __meta_Way = Ref{ProtoMeta}()
function meta(::Type{Way})
    ProtoBuf.metalock() do
        if !isassigned(__meta_Way)
            __meta_Way[] = target = ProtoMeta(Way)
            req = Symbol[:id]
            fnum = Int[1,2,3,4,8]
            pack = Symbol[:keys,:vals,:refs]
            wtype = Dict(:refs => :sint64)
            allflds = Pair{Symbol,Union{Type,String}}[:id => Int64, :keys => Base.Vector{UInt32}, :vals => Base.Vector{UInt32}, :info => Info, :refs => Base.Vector{Int64}]
            meta(target, Way, allflds, req, fnum, ProtoBuf.DEF_VAL, pack, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_Way[]
    end
end
function Base.getproperty(obj::Way, name::Symbol)
    if name === :id
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :keys
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :vals
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :info
        return (obj.__protobuf_jl_internal_values[name])::Info
    elseif name === :refs
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    else
        getfield(obj, name)
    end
end

const Relation_MemberType = (;[
    Symbol("NODE") => Int32(0),
    Symbol("WAY") => Int32(1),
    Symbol("RELATION") => Int32(2),
]...)

mutable struct Relation <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function Relation(; kwargs...)
        obj = new(meta(Relation), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct Relation
const __meta_Relation = Ref{ProtoMeta}()
function meta(::Type{Relation})
    ProtoBuf.metalock() do
        if !isassigned(__meta_Relation)
            __meta_Relation[] = target = ProtoMeta(Relation)
            req = Symbol[:id]
            fnum = Int[1,2,3,4,8,9,10]
            pack = Symbol[:keys,:vals,:roles_sid,:memids,:types]
            wtype = Dict(:memids => :sint64)
            allflds = Pair{Symbol,Union{Type,String}}[:id => Int64, :keys => Base.Vector{UInt32}, :vals => Base.Vector{UInt32}, :info => Info, :roles_sid => Base.Vector{Int32}, :memids => Base.Vector{Int64}, :types => Base.Vector{Int32}]
            meta(target, Relation, allflds, req, fnum, ProtoBuf.DEF_VAL, pack, wtype, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_Relation[]
    end
end
function Base.getproperty(obj::Relation, name::Symbol)
    if name === :id
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :keys
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :vals
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{UInt32}
    elseif name === :info
        return (obj.__protobuf_jl_internal_values[name])::Info
    elseif name === :roles_sid
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    elseif name === :memids
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int64}
    elseif name === :types
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Int32}
    else
        getfield(obj, name)
    end
end

mutable struct PrimitiveGroup <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function PrimitiveGroup(; kwargs...)
        obj = new(meta(PrimitiveGroup), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct PrimitiveGroup
const __meta_PrimitiveGroup = Ref{ProtoMeta}()
function meta(::Type{PrimitiveGroup})
    ProtoBuf.metalock() do
        if !isassigned(__meta_PrimitiveGroup)
            __meta_PrimitiveGroup[] = target = ProtoMeta(PrimitiveGroup)
            allflds = Pair{Symbol,Union{Type,String}}[:nodes => Base.Vector{Node}, :dense => DenseNodes, :ways => Base.Vector{Way}, :relations => Base.Vector{Relation}, :changesets => Base.Vector{ChangeSet}]
            meta(target, PrimitiveGroup, allflds, ProtoBuf.DEF_REQ, ProtoBuf.DEF_FNUM, ProtoBuf.DEF_VAL, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_PrimitiveGroup[]
    end
end
function Base.getproperty(obj::PrimitiveGroup, name::Symbol)
    if name === :nodes
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Node}
    elseif name === :dense
        return (obj.__protobuf_jl_internal_values[name])::DenseNodes
    elseif name === :ways
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Way}
    elseif name === :relations
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{Relation}
    elseif name === :changesets
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{ChangeSet}
    else
        getfield(obj, name)
    end
end

mutable struct PrimitiveBlock <: ProtoType
    __protobuf_jl_internal_meta::ProtoMeta
    __protobuf_jl_internal_values::Dict{Symbol,Any}
    __protobuf_jl_internal_defaultset::Set{Symbol}

    function PrimitiveBlock(; kwargs...)
        obj = new(meta(PrimitiveBlock), Dict{Symbol,Any}(), Set{Symbol}())
        values = obj.__protobuf_jl_internal_values
        symdict = obj.__protobuf_jl_internal_meta.symdict
        for nv in kwargs
            fldname, fldval = nv
            fldtype = symdict[fldname].jtyp
            (fldname in keys(symdict)) || error(string(typeof(obj), " has no field with name ", fldname))
            values[fldname] = isa(fldval, fldtype) ? fldval : convert(fldtype, fldval)
        end
        obj
    end
end # mutable struct PrimitiveBlock
const __meta_PrimitiveBlock = Ref{ProtoMeta}()
function meta(::Type{PrimitiveBlock})
    ProtoBuf.metalock() do
        if !isassigned(__meta_PrimitiveBlock)
            __meta_PrimitiveBlock[] = target = ProtoMeta(PrimitiveBlock)
            req = Symbol[:stringtable]
            val = Dict{Symbol,Any}(:granularity => 100, :lat_offset => 0, :lon_offset => 0, :date_granularity => 1000)
            fnum = Int[1,2,17,19,20,18]
            allflds = Pair{Symbol,Union{Type,String}}[:stringtable => StringTable, :primitivegroup => Base.Vector{PrimitiveGroup}, :granularity => Int32, :lat_offset => Int64, :lon_offset => Int64, :date_granularity => Int32]
            meta(target, PrimitiveBlock, allflds, req, fnum, val, ProtoBuf.DEF_PACK, ProtoBuf.DEF_WTYPES, ProtoBuf.DEF_ONEOFS, ProtoBuf.DEF_ONEOF_NAMES)
        end
        __meta_PrimitiveBlock[]
    end
end
function Base.getproperty(obj::PrimitiveBlock, name::Symbol)
    if name === :stringtable
        return (obj.__protobuf_jl_internal_values[name])::StringTable
    elseif name === :primitivegroup
        return (obj.__protobuf_jl_internal_values[name])::Base.Vector{PrimitiveGroup}
    elseif name === :granularity
        return (obj.__protobuf_jl_internal_values[name])::Int32
    elseif name === :lat_offset
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :lon_offset
        return (obj.__protobuf_jl_internal_values[name])::Int64
    elseif name === :date_granularity
        return (obj.__protobuf_jl_internal_values[name])::Int32
    else
        getfield(obj, name)
    end
end

export HeaderBlock, HeaderBBox, PrimitiveBlock, PrimitiveGroup, StringTable, Info, DenseInfo, ChangeSet, Node, DenseNodes, Way, Relation_MemberType, Relation


text:
struct Triangular{B <: NTuple{2, AbstractBoundary}} <: AbstractLattice
    dims::NTuple{2, Int}
    bcs::B
    translation_vectors::Dict{Int, Vector{Coordinate{2, Int}}}
    function Triangular{B}(dims::NTuple{2, Int}, bcs::B) where B <: NTuple{2, AbstractBoundary}
        check_boundaries(bcs)
        new{B}(dims, bcs, Dict{Int, Vector{Coordinate{2, Int}}}())
    end
end

for L in [:Honeycomb, :Kagome]

    @eval struct $L{B <: NTuple{2, AbstractBoundary}} <: AbstractLattice
        dims::NTuple{2, Int}
        bcs::B
        translation_vectors::Dict{Tuple{Int, Int}, Vector{Coordinate{3, Int}}}
        function $L{B}(dims::NTuple{2, Int}, bcs::B) where B <: NTuple{2, AbstractBoundary}
            check_boundaries(bcs)
            new{B}(dims, bcs, Dict{Tuple{Int, Int}, Vector{Coordinate{3, Int}}}())
        end
    end

    @eval hasunitcell(::$L) = Val{true}()
end

for L in [:Triangular, :Honeycomb, :Kagome]

    @eval ==(a::$L{B}, b::$L{B}) where B <: NTuple{2, AbstractBoundary} = (
        a.dims === b.dims
        && a.bcs === b.bcs
        && a.translation_vectors == b.translation_vectors
    )
    @eval ==(a::$L, b::$L) = false

    @eval ndims(::$L) = 2

    @eval metric(::$L) = Symmetric([1.0 0.5; 0.5 1.0])
    @eval ismetricdiag(::$L) = false

    @eval $L(dims::NTuple{2, Int}, bcs::NTuple{2, AbstractBoundary}) = $L{typeof(bcs)}(dims, bcs)
    @eval $L(dims::NTuple{2, Int}, bc::AbstractBoundary=Periodic()) = $L(dims, ntuple(x->bc, 2))

end

#############################
## Unitcell Properties
#############################

unitcellsize(::Honeycomb) = 2
unitcellsize(::Kagome) = 3

unitcell_translation_vectors(::Honeycomb) = [
    0.0 1/3;
    0.0 1/3
]

unitcell_translation_vectors(::Kagome) = [
    0.0 0.5 0.0;
    0.0 0.0 0.5
]

ext: jl


text:
Nx=200
Lx=1.0
u=1.0
dt=Lx/u/Nx/10
t_end=0.5*Lx/u
m=createMesh1D(Nx, Lx)
u_face=createFaceVariable(m, [u])
c_init=createCellVariable(m, 0.0)
c_init_up=createCellVariable(m, 0.0)
c=createCellVariable(m, 0.0)
c_upwind=createCellVariable(m, 0.0)

BC=createBC(m)
BC.left.a[:]=0.0
BC.left.b[:]=1.0
BC.left.c[:]=1.0
Mbc, RHSbc=boundaryConditionTerm(BC)

Mconv_up=convectionUpwindTerm(u_face)

flux_limiter=fluxLimiter("SUPERBEE")
for t in dt:dt:t_end
  Mt, RHSt=transientTerm(c_init, dt)
  Mt_up, RHSt_up=transientTerm(c_init_up, dt)
  c_upwind=solveLinearPDE(m, Mt_up+Mconv_up+Mbc, RHSt_up+RHSbc)
  for j in 1:10
    RHSconv=convectionTvdRHS(u_face, c, flux_limiter)
    c=solveLinearPDE(m, Mt+Mconv_up+Mbc, RHSt+RHSbc+RHSconv)
  end
  c_init=copyCell(c)
  c_init_up=copyCell(c_upwind)
end
Plots.plot(c.value[2:end-1])
Plots.plot!(c_upwind.value[2:end-1])


text:
using Documenter

makedocs(
    # modules = [TimesDates],
    sitename = "TimesDates",
    pages  = Any[
        "Overview"                     => "index.md",
        "Setup"                        => "setup.md",
        "Types"                        => "timedate_zone.md",
        "Examples"                     => "examples.md",
        "Design Notes"                 => "designnotes.md",
        "Acknowledgements"             => "acknowledgements.md"
        ]
    )

deploydocs(
    repo = "github.com/JeffreySarnoff/TimesDates.jl.git",
    target = "build"
)


text:
using Documenter, Lenny

makedocs(
    modules = [Lenny],
    format = :html,
    checkdocs = :exports,
    sitename = "Lenny.jl",
    pages = Any["index.md"]
)

deploydocs(
    repo = "github.com/dawbarton/Lenny.jl.git",
)


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

using Test
using Base.Meta
using Core: ReturnNode

"""
Helper to walk the AST and call a function on every node.
"""
function walk(func, expr)
    func(expr)
    if isa(expr, Expr)
        func(expr.head)
        for o in expr.args
            walk(func, o)
        end
    end
end

"""
Helper to test that every slot is in range after inlining.
"""
function test_inlined_symbols(func, argtypes)
    src, rettype = code_typed(func, argtypes)[1]
    nl = length(src.slotnames)
    ast = Expr(:block)
    ast.args = src.code
    walk(ast) do e
        if isa(e, Core.Slot)
            @test 1 <= e.id <= nl
        end
        if isa(e, Core.NewvarNode)
            @test 1 <= e.slot.id <= nl
        end
    end
end

# Test case 1:
# Make sure that all symbols are properly escaped after inlining
# https://github.com/JuliaLang/julia/issues/12620
@inline function test_inner(count)
    x = 1
    i = 0
    while i <= count
        y = x
        x = x + y
        i += 1
    end
end
function test_outer(a)
    test_inner(a)
end
test_inlined_symbols(test_outer, Tuple{Int64})

# Test case 2:
# Make sure that an error is thrown for the undeclared
# y in the else branch.
# https://github.com/JuliaLang/julia/issues/12620
@inline function foo_inl(x)
    if x
        y = 2
    else
        return y
    end
end
function bar12620()
    for i = 1:3
        foo_inl(i==1)
    end
end
@test_throws UndefVarError(:y) bar12620()

# issue #16165
@inline f16165(x) = (x = UInt(x) + 1)
g16165(x) = f16165(x)
@test g16165(1) === (UInt(1) + 1)

# issue #18948
f18948() = (local bar::Int64; bar=1.5)
g18948() = (local bar::Int32; bar=0x80000000)
@test_throws InexactError f18948()
@test_throws InexactError g18948()

# issue #21074
struct s21074
    x::Tuple{Int, Int}
end
@inline Base.getindex(v::s21074, i::Integer) = v.x[i]
@eval f21074() = $(s21074((1,2))).x[1]
let (src, _) = code_typed(f21074, ())[1]
    @test src.code[end] == ReturnNode(1)
end
@eval g21074() = $(s21074((1,2)))[1]
let (src, _) = code_typed(g21074, ())[1]
    @test src.code[end] == ReturnNode(1)
end

# issue #21311
counter21311 = Ref(0)
@noinline function update21311!(x)
    counter21311[] += 1
    x[] = counter21311[]
    return x
end
@noinline map21311(t::Tuple{Any}) = (update21311!(t[1]),)
@inline map21311(t::Tuple) = (update21311!(t[1]), map21311(Base.tail(t))...)
function read21311()
    xs = Ref(1), Ref(1)
    map21311(xs)
    return xs[1]
end
let a = read21311()
    @test a[] == 1
end

# issue #29083
f29083(;μ,σ) = μ + σ*randn()
g29083() = f29083(μ=2.0,σ=0.1)
let c = code_typed(g29083, ())[1][1].code
    # make sure no call to kwfunc remains
    @test !any(e->(isa(e,Expr) && ((e.head === :invoke && e.args[1].def.name === :kwfunc) ||
                                   (e.head === :foreigncall && e.args[1] === QuoteNode(:jl_get_keyword_sorter)))),
               c)
end

@testset "issue #19122: [no]inline of short func. def. with return type annotation" begin
    exf19122 = @macroexpand(@inline f19122()::Bool = true)
    exg19122 = @macroexpand(@noinline g19122()::Bool = true)
    @test exf19122.args[2].args[1].args[1] == :inline
    @test exg19122.args[2].args[1].args[1] == :noinline

    @inline f19122()::Bool = true
    @noinline g19122()::Bool = true
    @test f19122()
    @test g19122()
end

@testset "issue #27403: getindex is inlined with Union{Int,Missing}" begin
    function sum27403(X::AbstractArray)
        s = zero(eltype(X)) + zero(eltype(X))
        for x in X
            if !ismissing(x)
                s += x
            end
        end
        s
    end

    (src, _) = code_typed(sum27403, Tuple{Vector{Int}})[1]
    @test !any(x -> x isa Expr && x.head === :invoke, src.code)
end

# check that type.mutable can be fully eliminated
f_mutable_nothrow(s::String) = Val{typeof(s).mutable}
@test length(code_typed(f_mutable_nothrow, (String,))[1][1].code) == 1

# check that ifelse can be fully eliminated
function f_ifelse(x)
    a = ifelse(true, false, true)
    b = ifelse(a, true, false)
    return b ? x + 1 : x
end
# 2 for now because the compiler leaves a GotoNode around
@test_broken length(code_typed(f_ifelse, (String,))[1][1].code) <= 2

# Test that inlining of _apply properly hits the inference cache
@noinline cprop_inline_foo1() = (1, 1)
@noinline cprop_inline_foo2() = (2, 2)
function cprop_inline_bar(x...)
    if x === (1, 1, 1, 1)
        return x
    else
        # What you put here doesn't really matter,
        # the point is to prevent inlining when
        # x is not known to be (1, 1, 1, 1)
        println(stdout, "Hello")
        println(stdout, "World")
        println(stdout, "Hello")
        println(stdout, "World")
        println(stdout, "Hello")
        println(stdout, "World")
        println(stdout, "Hello")
        println(stdout, "World")
        println(stdout, "Hello")
        println(stdout, "World")
        println(stdout, "Hello")
        println(stdout, "World")
        return x
    end
    x
end

function cprop_inline_baz1()
    return cprop_inline_bar(cprop_inline_foo1()..., cprop_inline_foo1()...)
end
@test length(code_typed(cprop_inline_baz1, ())[1][1].code) == 1

function cprop_inline_baz2()
    return cprop_inline_bar(cprop_inline_foo2()..., cprop_inline_foo2()...)
end
@test length(code_typed(cprop_inline_baz2, ())[1][1].code) == 2

# Check that apply_type/TypeVar can be fully eliminated
function f_apply_typevar(T)
    NTuple{N, T} where N
    return T
end
@test length(code_typed(f_apply_typevar, (Type{Any},))[1][1].code) == 1

# check that div can be fully eliminated
function f_div(x)
	div(x, 1)
	return x
end
@test length(code_typed(f_div, (Int,))[1][1].code) == 1
# ...unless we div by an unknown amount
function f_div(x, y)
    div(x, y)
    return x
end
@test length(code_typed(f_div, (Int, Int))[1][1].code) > 1

f_identity_splat(t) = (t...,)
@test length(code_typed(f_identity_splat, (Tuple{Int,Int},))[1][1].code) == 1

# splatting one tuple into (,) plus zero or more empties should reduce
# this pattern appears for example in `fill_to_length`
f_splat_with_empties(t) = (()..., t..., ()..., ()...)
@test length(code_typed(f_splat_with_empties, (NTuple{200,UInt8},))[1][1].code) == 1

# check that <: can be fully eliminated
struct SomeArbitraryStruct; end
function f_subtype()
    T = SomeArbitraryStruct
    T <: Bool
end
let code = code_typed(f_subtype, Tuple{})[1][1].code
    @test length(code) == 1
    @test code[1] == ReturnNode(false)
end

# check that pointerref gets deleted if unused
f_pointerref(T::Type{S}) where S = Val(length(T.parameters))
let code = code_typed(f_pointerref, Tuple{Type{Int}})[1][1].code
    any_ptrref = false
    for i = 1:length(code)
        stmt = code[i]
        isa(stmt, Expr) || continue
        stmt.head === :call || continue
        arg1 = stmt.args[1]
        if arg1 === Base.pointerref || (isa(arg1, GlobalRef) && arg1.name === :pointerref)
            any_ptrref = true
        end
    end
    @test !any_ptrref
end

# Test that inlining can inline _applys of builtins/_applys on SimpleVectors
function foo_apply_apply_type_svec()
    A = (Tuple, Float32)
    B = Tuple{Float32, Float32}
    Core.apply_type(A..., B.types...)
end
let ci = code_typed(foo_apply_apply_type_svec, Tuple{})[1].first
    @test length(ci.code) == 1 && ci.code[1] == ReturnNode(NTuple{3, Float32})
end

# The that inlining doesn't drop ambiguity errors (#30118)
c30118(::Tuple{Ref{<:Type}, Vararg}) = nothing
c30118(::Tuple{Ref, Ref}) = nothing
b30118(x...) = c30118(x)

@test_throws MethodError c30118((Base.RefValue{Type{Int64}}(), Ref(Int64)))
@test_throws MethodError b30118(Base.RefValue{Type{Int64}}(), Ref(Int64))

# Issue #34900
f34900(x::Int, y) = x
f34900(x, y::Int) = y
f34900(x::Int, y::Int) = invoke(f34900, Tuple{Int, Any}, x, y)
let ci = code_typed(f34900, Tuple{Int, Int})[1].first
    @test length(ci.code) == 1 && isa(ci.code[1], ReturnNode) &&
        ci.code[1].val.n == 2
end

@testset "check jl_ir_flag_inlineable for inline macro" begin
    @test ccall(:jl_ir_flag_inlineable, Bool, (Any,), first(methods(@inline x -> x)).source)
    @test !ccall(:jl_ir_flag_inlineable, Bool, (Any,), first(methods( x -> x)).source)
    @test ccall(:jl_ir_flag_inlineable, Bool, (Any,), first(methods(@inline function f(x) x end)).source)
    @test !ccall(:jl_ir_flag_inlineable, Bool, (Any,), first(methods(function f(x) x end)).source)
end

const _a_global_array = [1]
f_inline_global_getindex() = _a_global_array[1]
let ci = code_typed(f_inline_global_getindex, Tuple{})[1].first
    @test any(x->(isexpr(x, :call) && x.args[1] === GlobalRef(Base, :arrayref)), ci.code)
end

# Issue #29114 & #36087 - Inlining of non-tuple splats
f_29115(x) = (x...,)
@test @allocated(f_29115(1)) == 0
@test @allocated(f_29115(1=>2)) == 0
let ci = code_typed(f_29115, Tuple{Int64})[1].first
    @test length(ci.code) == 2 && isexpr(ci.code[1], :call) &&
        ci.code[1].args[1] === GlobalRef(Core, :tuple)
end
let ci = code_typed(f_29115, Tuple{Pair{Int64, Int64}})[1].first
    @test length(ci.code) == 4 && isexpr(ci.code[1], :call) &&
        ci.code[end-1].args[1] === GlobalRef(Core, :tuple)
end

# Issue #37182 & #37555 - Inlining of pending nodes
function f37555(x::Int; kwargs...)
    @assert x < 10
    +(x, kwargs...)
end
@test f37555(1) == 1

# Test that we can inline small constants even if they are not isbits
struct NonIsBitsDims
    dims::NTuple{N, Int} where N
end
NonIsBitsDims() = NonIsBitsDims(())
let ci = code_typed(NonIsBitsDims, Tuple{})[1].first
    @test length(ci.code) == 1 && isa(ci.code[1], ReturnNode) &&
        ci.code[1].val.value == NonIsBitsDims()
end

struct NonIsBitsDimsUndef
    dims::NTuple{N, Int} where N
    NonIsBitsDimsUndef() = new()
end
@test Core.Compiler.is_inlineable_constant(NonIsBitsDimsUndef())
@test !Core.Compiler.is_inlineable_constant((("a"^1000, "b"^1000), nothing))

max_stars_repo_name: tckmn/julia

max_issues_repo_name: tckmn/julia

max_forks_repo_name: tckmn/julia


text:
using MPI, Test

@testset "Atmos/Dycore" begin
  # The code below was modified from the MPI.jl file runtests.jl
  #
  # Code coverage command line options; must correspond to src/julia.h
  # and src/ui/repl.c
  JL_LOG_NONE = 0
  JL_LOG_USER = 1
  JL_LOG_ALL = 2
  coverage_opts = Dict{Int, String}(JL_LOG_NONE => "none",
                                    JL_LOG_USER => "user",
                                    JL_LOG_ALL => "all")
  coverage_opt = coverage_opts[Base.JLOptions().code_coverage]
  testdir = dirname(@__FILE__)

  for (n, f) in [(1, "tracers_test.jl")
                 (1, "rising_thermal_bubble.jl")
                ]
    cmd =  `mpiexec -n $n $(Base.julia_cmd()) --startup-file=no --project=$(Base.active_project()) --code-coverage=$coverage_opt $(joinpath(testdir, f))`
    @info "Running MPI test..." n f cmd
    # Running this way prevents:
    #   Atmos/Dycore | No tests
    # since external tests are not returned as passed/fail
    @test (run(cmd); true)
  end
end

max_stars_repo_name: jkozdon/CLIMA

max_issues_repo_name: jkozdon/CLIMA

max_forks_repo_name: jkozdon/CLIMA


text:
"""
    pathfinder(logp, ∇logp, θ₀::AbstractVector{<:Real}, ndraws::Int; kwargs...)

Find the best multivariate normal approximation encountered while optimizing `logp`.

The multivariate normal approximation returned is the one that maximizes the evidence lower
bound (ELBO), or equivalently, minimizes the KL divergence between

# Arguments
- `logp`: a callable that computes the log-density of the target distribution.
- `∇logp`: a callable that computes the gradient of `logp`.
- `θ₀`: initial point of length `dim` from which to begin optimization
- `ndraws`: number of approximate draws to return

# Keywords
- `rng::Random.AbstractRNG`: The random number generator to be used for drawing samples
- `optimizer::Optim.AbstractOptimizer`: Optimizer to be used for constructing trajectory.
Defaults to `Optim.LBFGS(; m=$DEFAULT_HISTORY_LENGTH, linesearch=$DEFAULT_LINE_SEARCH)`.
- `history_length::Int=$DEFAULT_HISTORY_LENGTH`: Size of the history used to approximate the
inverse Hessian. This should only be set when `optimizer` is not an `Optim.LBFGS`.
- `ndraws_elbo::Int=5`: Number of draws used to estimate the ELBO
- `kwargs...` : Remaining keywords are forwarded to `Optim.Options`.

# Returns
- `q::Distributions.MvNormal`: ELBO-maximizing multivariate normal distribution
- `ϕ::AbstractMatrix{<:Real}`: draws from multivariate normal with size `(dim, ndraws)`
- `logqϕ::Vector{<:Real}`: log-density of multivariate normal at columns of `ϕ`
"""
function pathfinder(
    logp,
    ∇logp,
    θ₀,
    ndraws;
    rng::Random.AbstractRNG=Random.default_rng(),
    optimizer=DEFAULT_OPTIMIZER,
    history_length::Int=optimizer isa Optim.LBFGS ? optimizer.m : DEFAULT_HISTORY_LENGTH,
    ndraws_elbo::Int=5,
    kwargs...,
)
    # compute trajectory
    θs, logpθs, ∇logpθs = maximize_with_trace(logp, ∇logp, θ₀, optimizer; kwargs...)
    L = length(θs) - 1
    @assert L + 1 == length(logpθs) == length(∇logpθs)
    # fit mv-normal distributions to trajectory
    qs = fit_mvnormals(θs, ∇logpθs; history_length=history_length)
    # find ELBO-maximizing distribution
    lopt, elbo, ϕ, logqϕ = maximize_elbo(rng, logp, qs[2:end], ndraws_elbo)
    @info "Optimized for $L iterations. Maximum ELBO of $(round(elbo; digits=2)) reached at iteration $lopt."

    # get parameters of ELBO-maximizing distribution
    q = qs[lopt + 1]

    # reuse existing draws; draw additional ones if necessary
    if ndraws_elbo < ndraws
        ϕnew, logqϕnew = rand_and_logpdf(rng, q, ndraws - ndraws_elbo)
        ϕ = hcat(ϕ, ϕnew)
        append!(logqϕ, logqϕnew)
    elseif ndraws < ndraws_elbo
        ϕ = ϕ[:, 1:ndraws]
        logqϕ = logqϕ[1:ndraws]
    end

    return q, ϕ, logqϕ
end


text:
println("Loading NGWP...")
using NGWP
using Test, JLD, MAT, Plots, LightGraphs, LinearAlgebra, SparseArrays, MTSG

#####################################################
# 1. Testing PC-NGWP and VM-NGWP on sunflower graph #
#####################################################
println("Testing NGWP on sunflower barbara eye signal")
## Build weighted sunflower graph for test
G, L, X = SunFlowerGraph(N = 400); N = nv(G)
𝛌, 𝚽 = eigen(Matrix(L)); standardize_eigenvectors!(𝚽)
Q = incidence_matrix(G; oriented = true)
W = 1.0 * adjacency_matrix(G)
edge_weight = [e.weight for e in edges(G)]

## Build Dual Graph by DAG metric
distDAG = eigDAG_Distance(𝚽, Q, N; edge_weight = edge_weight)
Gstar_Sig = dualgraph(distDAG)
G_Sig = GraphSig(W, xy = X)
GP_dual = partition_tree_fiedler(Gstar_Sig; swapRegion = false)
GP_primal = pairclustering(𝚽, GP_dual);

VM_NGWP = vm_ngwp(𝚽, GP_dual)
PC_NGWP = pc_ngwp(𝚽, GP_dual, GP_primal)

## use barbara eye signal for testing
using MAT
f = matread(joinpath(@__DIR__, "datasets",
                "sunflower_barbara_voronoi.mat"))["f_eye_voronoi"]
G_Sig.f = reshape(f, (N, 1))

DVEC = getall_expansioncoeffs(G_Sig, GP_dual, VM_NGWP, PC_NGWP, 𝚽)
DVEC_truth = JLD.load(joinpath(@__DIR__,
                        "datasets", "sunflower_barbara_feye_DVEC.jld"), "DVEC")

@testset "PC-NGWP Test" begin
    @test norm(DVEC[end-1] - DVEC_truth[end - 1]) / norm(DVEC_truth[end - 1]) < 1e-8
end

@testset "VM-NGWP Test" begin
    @test norm(DVEC[end] - DVEC_truth[end]) / norm(DVEC_truth[end]) < 1e-8
end


text:
using AbstractPattern
using AbstractPattern.BasicPatterns
using Test

const backend = MK(RedyFlavoured)

# function test()

# pairit(x) = :(complex_func1($x), complex_func2($x))

# code = backend(
#     :val,
#     [
#         or(
#             literal(1),
#             literal("string"),
#         ) => :a,

#         and(literal(2), P_capture(:a)) => :b,

#         guard((_, _, _) -> :some_cond) => :c,

#         and(
#             P_type_of(Symbol),
#             P_bind(:c, :((val, val)))
#         ) => :d,

#         literal(3) => :e,

#         P_vector([
#             literal(1),
#             and(P_type_of(AbstractFloat),  P_capture(:a)),
#             or(literal(1), P_type_of(Real))
#         ]) =>  :f,

#         P_tuple([
#             literal(5),
#             P_type_of(AbstractString),
#             P_tuple([
#                 P_type_of(Int),
#                 P_type_of(Int)
#             ])
#         ]) =>  :g,

#         # pairit(x) = :(complex_func1($x), complex_func2($x))
#         P_slow_view(pairit,
#             P_tuple([P_type_of(Int), P_type_of(Int)])
#         ) => :h,
#         P_slow_view(pairit,
#             P_tuple([P_type_of(String), P_type_of(Int)])
#         ) => :i
#     ]
# )

# end

# println(test())

untyped_view1(x) = :(f1($x))
untyped_view2(x) = :(f2($x))

code = backend(
    :val,
    [
        LineNumberNode(@__LINE__),
        and(
            P_slow_view(untyped_view2,
                [P_type_of(Int), P_type_of(Int)] |> P_tuple,
            ),
            P_slow_view(untyped_view1, literal(true)),

        ) => :h,
        P_slow_view(untyped_view2,
            [
                P_type_of(String),
                P_slow_view(untyped_view1, literal(false))
            ] |> P_tuple,
        ) => :i
    ]
)
println(code)

max_stars_repo_name: KristofferC/AbstractPattern

max_issues_repo_name: KristofferC/AbstractPattern

max_forks_repo_name: KristofferC/AbstractPattern


text:
# call stack:
# run/replay
#   > calls apply(GreensIterator, group, dqmc, model, sweep)
#       > calls measure(LatticeIterator, measurement, dqmc, model, sweep, GreensResults...)
#           > calls apply!(LatticeIterator, measurement, dqmc, model, GreensResults...)
#               > calls kernel(measurement, dqmc, model, MaskResult..., GreensResults...)

# Does:
# creates groups, runs simulation
#   > resolves GreensIterator (calclates G, Gkk, Gkl or whatever)
#       > resolves sweep (skipping, frequency based measurements) and commits
#           > resolves Mask (to Greens indices)
#               > calculate element from Wicks theorem

# TODO
# - deprecate old measurement types
# - _save/_load to autoconvert
# - re-add Greens -> Occupation maybe?

struct DQMCMeasurement{
        GreensIterator, LatticeIterator, F <: Function, OT, T
    } <: AbstractMeasurement
    kernel::F
    observable::OT
    output::T
end

function DQMCMeasurement{GI, LI}(kernel::FT, observable::OT, output::T) where {GI, LI, FT, OT, T}
    DQMCMeasurement{GI, LI, FT, OT, T}(kernel, observable, output)
end

function Measurement(
        dqmc, model, GreensIterator, LatticeIterator, kernel::F;
        capacity = _default_capacity(dqmc), eltype = geltype(dqmc),
        shape = _get_shape(dqmc, model, LatticeIterator)
    ) where F
    temp = shape === nothing ? zero(eltype) : zeros(eltype, shape)
    obs = LogBinner(temp, capacity=capacity)
    DQMCMeasurement{GreensIterator, LatticeIterator}(kernel, obs, temp)
end

################################################################################
### DQMCMeasurement utilities
################################################################################

function Base.show(io::IO, ::MIME"text/plain", m::DQMCMeasurement{GI, LI}) where {GI, LI}
    max = capacity(m.observable)
    current = length(m.observable)
    print(io, "[$current/$max] DQMCMeasurement{$GI, $LI}($(m.kernel))")
end

MonteCarloObservable.mean(m::DQMCMeasurement) = mean(m.observable)
MonteCarloObservable.var(m::DQMCMeasurement) = var(m.observable)
MonteCarloObservable.std_error(m::DQMCMeasurement) = std_error(m.observable)
MonteCarloObservable.tau(m::DQMCMeasurement) = tau(m.observable)
Base.length(m::DQMCMeasurement) = length(m.observable)
Base.isempty(m::DQMCMeasurement) = isempty(m.observable)
Base.empty!(m::DQMCMeasurement) = empty!(m.observable)

"""
    _default_capacity(mc::DQMC)

Returns 100_000 if the configuration measurement is not available or empty, and
the length of the configuration measurement if it is not.

This should be useful for `run` - `setup measurements` - `replay` workflow,
returning the exact capacity needed.
"""
function _default_capacity(mc::DQMC)
    k = if isdefined(mc, :measurements)
        findfirst(v -> v isa ConfigurationMeasurement, mc.measurements)
    else
        nothing
    end
    if k === nothing
        return 100_000
    else
        N = length(mc.measurements[k].obs)
        return N == 0 ? 100_000 : N
    end
end

# _get_shape(model, ::Nothing) = (length(lattice(model)),)
_get_shape(model, mask::RawMask) = (mask.nsites, mask.nsites)
_get_shape(model, mask::DistanceMask) = length(mask)

_get_shape(mc, model, LI::Type) = _get_shape(model, LI(mc, model))
_get_shape(model, ::Nothing) = nothing
_get_shape(model, ::EachSite) = length(lattice(model))
_get_shape(model, ::EachSiteAndFlavor) = nflavors(model) * length(lattice(model))
_get_shape(model, ::EachSitePair) = (length(lattice(model)), length(lattice(model)))
_get_shape(model, iter::EachSitePairByDistance) = ndirections(iter)
_get_shape(model, iter::EachLocalQuadByDistance) = ndirections(iter)
_get_shape(model, iter::EachLocalQuadBySyncedDistance) = ndirections(iter)

# Some type piracy to make things easier
Base.Nothing(::DQMC, ::Model) = nothing

# To identify the requirement of equal-time Greens functions
struct Greens <: AbstractGreensIterator end
Greens(::DQMC, ::Model)= Greens()

struct GreensAt{k, l} <: AbstractUnequalTimeGreensIterator end
GreensAt(l::Integer) = GreensAt{l, l}
GreensAt(k::Integer, l::Integer) = GreensAt{k, l}
GreensAt{k, l}(::DQMC, ::Model) where {k, l} = GreensAt{k, l}()

# maybe we want to some stuff to get custom indexing?
# How would this integrate with greens iterators? maybe:
# greens(mc::DQMC, model::Model) = greens(mc)

requires(::AbstractMeasurement) = (Nothing, Nothing)
requires(::DQMCMeasurement{GI, LI}) where {GI, LI} = (GI, LI)

@bm function generate_groups(mc, model, measurements)
    # maybe instead:
    requirements = requires.(measurements)
    GIs = tuple(unique(first.(requirements))...)
    LIs = tuple(unique(last.(requirements))...)
    lattice_iterators = map(T -> T(mc, model), LIs)

    if any(T -> T <: AbstractUnequalTimeGreensIterator, GIs)
        initialize_stack(mc, mc.ut_stack)
    end

    # greens_iterator => [
    #     (lattice_iterator, measurement),
    #     (lattice_iterator, measurement),
    #     ...
    # ]
    map(enumerate(GIs)) do (i, G)
        ms = filter(m -> G == requires(m)[1], measurements)
        group = map(ms) do m
            LI = requires(m)[2]
            j = findfirst(==(LI), LIs)
            @assert j !== nothing
            (lattice_iterators[j], m)
        end
        G(mc, model) => group
    end
end

lattice_iterator(::DQMCMeasurement{GI, LI}, mc, model) where {GI, LI} = LI(mc, model)

# TODO
# Saving the kernel function as a symbol is kinda risky because the function
# definition is not guaranteed to be the same
# though maybe that's a good thing - changing function definitions doesn't
# break stuff this way
function _save(file::JLDFile, m::DQMCMeasurement{GI, LI}, key::String) where {GI, LI}
    write(file, "$key/VERSION", 1)
    write(file, "$key/type", DQMCMeasurement)
    write(file, "$key/GI", GI)
    write(file, "$key/LI", LI)
    # maybe add module for eval?
    write(file, "$key/kernel", Symbol(m.kernel))
    write(file, "$key/obs", m.observable)
    write(file, "$key/output", m.output)
end

# TODO
# I think eval will fail if kernel is not defined in MonteCarlo
function _load(data, ::Type{T}) where {T <: DQMCMeasurement}
    kernel = try
        eval(data["kernel"])
    catch e
        @warn "Failed to load kernel in module MonteCarlo." exception=e
        missing_kernel
    end
    DQMCMeasurement{data["GI"], data["LI"]}(kernel, data["obs"], data["output"])
end

missing_kernel(args...) = error("kernel couldn't be loaded.")

################################################################################
### Greens function related
################################################################################

@bm function apply!(::Nothing, combined::Vector{<: Tuple}, mc::DQMC, model, sweep)
    for (lattice_iterator, measurement) in combined
        # Clear output if necessary
        prepare!(lattice_iterator, model, measurement)
        # Write measurement to ouput
        measure!(lattice_iterator, measurement, mc, model, sweep)
        # Finalize output and commit
        finish!(lattice_iterator, model, measurement)
    end

    nothing
end

@bm function apply!(::Greens, combined::Vector{<: Tuple}, mc::DQMC, model, sweep)
    G = greens!(mc)
    for (lattice_iterator, measurement) in combined
        prepare!(lattice_iterator, model, measurement)
        measure!(lattice_iterator, measurement, mc, model, sweep, G)
        finish!(lattice_iterator, model, measurement)
    end
    nothing
end

@bm function apply!(::GreensAt{k, l}, combined::Vector{<: Tuple}, mc::DQMC, model, sweep) where {k, l}
    G = greens!(mc, k, l)
    for (lattice_iterator, measurement) in combined
        prepare!(lattice_iterator, model, measurement)
        measure!(lattice_iterator, measurement, mc, model, sweep, G)
        finish!(lattice_iterator, model, measurement)
    end
    nothing
end

@bm function apply!(iter::CombinedGreensIterator, combined::Vector{<: Tuple}, mc::DQMC, model, sweep)
    for (lattice_iterator, measurement) in combined
        prepare!(lattice_iterator, model, measurement)
    end

    G00 = greens!(mc)
    for (G0l, Gl0, Gll) in iter
        for (lattice_iterator, measurement) in combined
            measure!(lattice_iterator, measurement, mc, model, sweep, (G00, G0l, Gl0, Gll))
        end
    end

    for (lattice_iterator, measurement) in combined
        finish!(lattice_iterator, model, measurement, mc.p.delta_tau)
    end
    nothing
end

################################################################################
### measure!
################################################################################

@bm function measure!(lattice_iterator, measurement, mc::DQMC, model, sweep, packed_greens)
    # ignore sweep
    apply!(lattice_iterator, measurement, mc, model, packed_greens)
    nothing
end

# Lattice irrelevant
@bm function measure!(::Nothing, measurement, mc::DQMC, model, sweep, packed_greens)
    push!(measurement.observable, measurement.kernel(mc, model, packed_greens))
    nothing
end

################################################################################
### LatticeIterator related
################################################################################

# If LatticeIterator is Nothing, then things should be handled in measure!
prepare!(::Nothing, model, m) = nothing
prepare!(::AbstractLatticeIterator, model, m) = m.output .= zero(eltype(m.output))

finish!(::Nothing, model, m) = nothing # handled in measure!
finish!(::AbstractLatticeIterator, model, m) = push!(m.observable, m.output)
function finish!(::AbstractLatticeIterator, model, m, factor)
    m.output .*= factor
    push!(m.observable, m.output)
end
function finish!(::EachSitePairByDistance, model, m, factor=1.0)
    m.output .*= factor / length(lattice(model))
    push!(m.observable, m.output)
end
function finish!(::EachLocalQuadByDistance, model, m, factor=1.0)
    m.output .*= factor / length(lattice(model))
    push!(m.observable, m.output)
end
function finish!(::EachLocalQuadBySyncedDistance, model, m, factor=1.0)
    m.output .*= factor / length(lattice(model))
    push!(m.observable, m.output)
end

# Call kernel for each site (linear index)
@bm function apply!(iter::EachSiteAndFlavor, measurement, mc::DQMC, model, packed_greens)
    for i in iter
        measurement.output[i] += measurement.kernel(mc, model, i, packed_greens)
    end
    nothing
end

# Call kernel for each site (linear index)
@bm function apply!(iter::EachSite, measurement, mc::DQMC, model, packed_greens)
    for i in iter
        measurement.output[i] += measurement.kernel(mc, model, i, packed_greens)
    end
    nothing
end

# Call kernel for each pair (src, trg) (Nsties² total)
@bm function apply!(iter::EachSitePair, measurement, mc::DQMC, model, packed_greens)
    for (i, j) in iter
        measurement.output[i, j] += measurement.kernel(mc, model, i, j, packed_greens)
    end
    nothing
end

# Call kernel for each pair (site, site) (i.e. on-site)
@bm function apply!(iter::OnSite, measurement, mc::DQMC, model, packed_greens)
    for (i, j) in iter
        measurement.output[i] += measurement.kernel(mc, model, i, j, packed_greens)
    end
    nothing
end

# Call kernel for each pair (src, trg) and sum those that point in the same direction
@bm function apply!(iter::EachSitePairByDistance, measurement, mc::DQMC, model, packed_greens)
    @inbounds for (dir, src, trg) in iter
        measurement.output[dir] += measurement.kernel(mc, model, src, trg, packed_greens)
    end
    nothing
end

# Call kernel for each pair (src1, trg1, src2, trg) and sum those that have the
# same `dir12 = pos[src2] - pos[src1]`, `dir1 = pos[trg1] - pos[src1]` and
# `dir2 = pos[trg2] - pos[src2]`
@bm function apply!(iter::EachLocalQuadByDistance, measurement, mc::DQMC, model, packed_greens)
    # lin is a linear index for (dir12, dir1, dir2)
    for (lin, src1, trg1, src2, trg2) in iter
        measurement.output[lin] += measurement.kernel(
            mc, model, src1, trg1, src2, trg2, packed_greens
        )
    end
    nothing
end

# Call kernel for each pair (src1, trg1, src2, trg) and sum those that have the
# same `dir12 = pos[src2] - pos[src1]` and
# `dir1 = pos[trg1] - pos[src1] = dir2 = pos[trg2] - pos[src2] = dir_ii`
@bm function apply!(iter::EachLocalQuadBySyncedDistance, measurement, mc::DQMC, model, packed_greens)
    # lin is a linear index for (dir12, dir_ii)
    for (lin, src1, trg1, src2, trg2) in iter
        measurement.output[lin] += measurement.kernel(
            mc, model, src1, trg1, src2, trg2, packed_greens
        )
    end
    nothing
end

include("measurements.jl")
include("deprecated.jl")


text:
"""
    BNKeyType

Keys can hold value of type T or nothing.
"""
BNKeyType = Union{T, Nothing} where T

"""
    BNode

Binary node.
"""
mutable struct BNode{T}
    key::BNKeyType{T}
    left::BNode{T}
    right::BNode{T}
    BNode{T}() where T = (x = new(); x)
    BNode{T}(k::BNKeyType{T}) where T = (x = new(); x.key = k; x)
    BNode{T}(k::BNKeyType{T}, l::BNode{T}, r::BNode{T}) where T = new(k, l, r)
end

"""
    RMHeap

Randomised meldable heaps type.
"""
mutable struct RMHeap{T}
    root::BNode{T}
    N::Int
    RMHeap{T}() where T = new(BNode{T}(), 0)
end

"""
    add!(rmheap, x)

Add to the heap.
"""
function add!(rmheap::RMHeap{T}, x::T) where T
    rmheap.root = meld(rmheap.root, BNode{T}(x))
    rmheap.N += 1
    return nothing
end

"""
    remove(rmheap)

Remove minimum element from the heap.
"""
function remove!(rmheap::RMHeap{T}) where T
    x = rmheap.root.key
    rmheap.N -= 1
    if !isdefined(rmheap.root, :left) && !isdefined(rmheap.root, :right)
        rmheap.root = BNode{T}()
        return x
    end
    if !isdefined(rmheap.root, :left)
        rmheap.root = rmheap.root.right
        return x
    end
    if !isdefined(rmheap.root, :right)
        rmheap.root = rmheap.root.left
        return x
    end
    rmheap.root = meld(rmheap.root.left, rmheap.root.right)

    return x
end

"""
    meld!(rmheap, other)

Meld two heaps.
"""
function meld!(rmheap::RMHeap{T}, other::RMHeap{T}) where T
    rmheap.root = meld(rmheap.root, other.root)
    other.root = BNode{T}()
    return nothing
end

"""
    meld(node1, node2)

Meld two nodes.
"""
function meld(node1::BNode{T}, node2::BNode{T}) where T
    if node1.key === nothing return node2 end
    if node2.key === nothing return node1 end

    if node1.key > node2.key
        tmp = node1
        node1 = node2
        node2 = tmp
    end

    if rand((0, 1)) == 1
        if !isdefined(node1, :left) node1.left = BNode{T}() end
        node1.left = meld(node1.left, node2)
    else
        if !isdefined(node1, :right) node1.right = BNode{T}() end
        node1.right = meld(node1.right, node2)
    end

    return node1
end

"""
    getmin(rmheap)

Get the minimum element, don't remove.
"""
getmin(rmheap::RMHeap) = rmheap.root.key

"""
    length(rmheap)
    size(rmheap)

Get size of priority queue.
"""
length(rmheap::RMHeap) = rmheap.N
size(rmheap::RMHeap) = length(rmheap)

# iterate
"""
    eltype(type)

Type of elements generated by iteration over a linked list.
"""
eltype(::RMHeap{T}) where T = T

"""
    isempty(rmheap)

Check if priority queue is empty.
"""
isempty(rmheap::RMHeap) = rmheap.N == 0

"""
    MinPQ

Binary heap minimum priority queue.
"""
mutable struct MinPQ{T}
    root::Vector{T}
    N::Int
    MinPQ{T}(n::Int) where T = new(Vector{T}(undef, n), n)
    MinPQ{T}() where T = new(Vector{T}(undef, 0), 0)
end

"""
    add(pq, x)

Add element to the priority queue.
"""
function add!(pq::MinPQ{T}, x::T) where T
    pq.N += 1
    push!(pq.root, x)
    swim!(pq, pq.N)

    return nothing
end

"""
    remove!(pq)

Remove and get the smallst element from the priority queue.
"""
function remove!(pq::MinPQ)
    if pq.N == 0 return nothing end

    min = pq.root[1]
    pq.root[1] = pq.root[pq.N]
    pq.root[pq.N] = min
    pop!(pq.root)
    pq.N -= 1
    sink!(pq, 1)

    return min
end

"""
    swim!(pq)

Increase an elements priority in the priority queue.
"""
function swim!(pq::MinPQ, k::Int)
    while k > 1 && greater(pq, k ÷ 2, k)
        tmp = pq.root[k ÷ 2]
        pq.root[k ÷ 2] = pq.root[k]
        pq.root[k] = tmp
        k = k ÷ 2
    end

    return nothing
end

"""
    sink!(pq, k)

Lower a node in the priority queue.
"""
function sink!(pq::MinPQ, k::Int)
    while 2k <= pq.N
        j = 2k
        if j < pq.N && greater(pq, j, j + 1) j += 1 end
        if !greater(pq, k, j) break end
        tmp = pq.root[k]
        pq.root[k] = pq.root[j]
        pq.root[j] = tmp
        k = j
    end

    return nothing
end

"""
    isempty(pq)

Check if priority queue is empty.
"""
isempty(pq::MinPQ) = pq.N == 0

"""
    getmin(pq)

Get the smallest element from priority queue, don't remove.
"""
getmin(pq::MinPQ) = pq.root[1]

"""
    greater(pq, i, j)

Check if element i is larger than element j in priority queue.
"""
greater(pq::MinPQ, i::Int, j::Int) = pq.root[i] > pq.root[j]

"""
    length(pq)
    size(pq)

Get size of priority queue.
"""
length(pq::MinPQ) = pq.N
size(pq::MinPQ) = length(pq)

"""
    isheap(heap)
    isheap(node)
    isheap(pq)

Check that heap invariant holds.
"""
isheap(heap::RMHeap) = isheap(heap.root)
function isheap(node::BNode)
    if isdefined(node, :left)
        if node.key > node.left.key
            return false
        end
        isheap(node.left)
    end
    if isdefined(node, :right)
        if node.key > node.right.key
            return false
        end
        isheap(node.right)
    end

    return true
end

function isheap(pq::MinPQ)
    if isempty(pq) return true end
    data = pq.root

    i = 1
    while 2i <= pq.N
        if !(data[i] <= data[2i]) return false end
        if 2i + 1 > pq.N break end
        if !(data[i] <= data[2i + 1]) return false end
        i += 1
    end

    return true
end


text:
using CUDA, FixedEffectModels, DataFrames, CSV, Test

df = DataFrame(CSV.File(joinpath(dirname(pathof(FixedEffectModels)), "../dataset/Cigar.csv")))

##############################################################################
##
## Printing Results
##
##############################################################################

model = @formula Sales ~ NDI
result = reg(df, model)
show(result)
predict(result, df)
residuals(result, df)
@test responsename(result) == "Sales"

model = @formula Sales ~ CPI + (Price ~ Pimin)
result = reg(df, model)
coeftable(result)
show(result)

predict(result, df)
residuals(result, df)
@test  nobs(result) == 1380
@test vcov(result)[1] ≈ 3.5384578251636785

model = @formula Sales ~ Price + fe(State)
result = reg(df, model)
show(result)
model = @formula Sales ~ CPI + (Price ~ Pimin) + fe(State)
result = reg(df, model)
show(result)

## Tests for predict method
# Test that only DataFrame is accepted
@test_throws "Predict requires an input of type DataFrame" predict(result, Matrix(df))

# Test that predicting from model without saved FE test throws
@test_throws "No estimates for fixed effects found. Model needs to be estimated with save = :fe or :all for prediction to work." predict(result, df)

# Test basic functionality - adding 1 to price should increase prediction by coef
model = @formula Sales ~ Price + fe(State)
result = reg(df, model, save = :fe)
x = predict(result, DataFrame(Price = [1.0, 2.0], State = [1, 1]))
@test last(x) - first(x) ≈ only(result.coef)

# Missing variables in covariates should yield missing prediction
x = predict(result, DataFrame(Price = [1.0, missing], State = [1, 1]))
@test ismissing(last(x))

# Missing variables in fixed effects should yield missing prediction
x = predict(result, DataFrame(Price = [1.0, 2.0], State = [1, missing]))
@test ismissing(last(x))

# Fixed effect levels not in the estimation data should yield missing prediction
x = predict(result, DataFrame(Price = [1.0, 2.0], State = [1, 111]))
@test ismissing(last(x))

##############################################################################
##
## Saved Residuals
##
##############################################################################

model = @formula Sales ~ Price
result = reg(df, model)
@test residuals(result, df)[1:10] ≈ [-39.2637, -37.48801, -34.38801, -36.09743, -36.97446, -43.15547, -41.22573, -40.83648, -34.52427, -28.91617] atol = 1e-4
@test r2(result) ≈ 0.0968815737054879 atol = 1e-4
@test adjr2(result) ≈ 0.0962261902321246 atol = 1e-4
@test result.nobs == 1380
@test result.F ≈ 147.8242550248069 atol= 1e-4

#weights
model = @formula Sales ~ CPI
result = reg(df, model, weights = :Pop)
@test residuals(result, df)[1:3] ≈ [ -35.641449,  -34.0611538,  -30.860784] atol = 1e-4

# iv
model = @formula Sales ~ CPI + (Price ~ Pimin)
result = reg(df, model)
@test residuals(result, df)[1:3] ≈ [ -33.047390, -30.9518422, -28.1371048] atol = 1e-4

# iv with exo after endo
model = @formula Sales ~ (Price ~ Pimin) + CPI
result = reg(df, model)
@test residuals(result, df)[1:3] ≈ [ -33.047390, -30.9518422, -28.1371048] atol = 1e-4

# iv and weights
model = @formula Sales ~ CPI + (Price ~ Pimin)
result = reg(df, model, weights = :Pop)
@test residuals(result, df)[1:3] ≈ [ -30.2284549, -28.09507, -25.313248] atol = 1e-4

# iv, weights and subset of states
model = @formula Sales ~ CPI + (Price ~ Pimin)
result = reg(df, model, subset = df.State .<= 30, weights = :Pop)
@test residuals(result, df)[1:3] ≈ [ -34.081720, -31.906020, -29.131738] atol = 1e-4

# fixed effects
model = @formula Sales ~ Price + fe(State)
result = reg(df, model, save = true)
@test residuals(result)[1:3] ≈ [-22.08499, -20.33318, -17.23318] atol = 1e-4
@test result.nobs == 1380
@test r2(result) ≈ 0.7682403747044817 atol = 1e-4
@test adjr2(result) ≈ 0.7602426682051615 atol = 1e-4
@test result.F ≈ 458.4582526109375 atol = 1e-4

# fixed effects and weights
model = @formula Sales ~ Price + fe(State)
result = reg(df, model,  weights = :Pop, save = true)
@test residuals(result)[1:3] ≈ [ -23.413793, -21.65289, -18.55289] atol = 1e-4

# fixed effects and iv
#TO CHECK WITH IVREGHDFE, NO SUPPORT RIGHT NOW
model = @formula Sales ~ CPI + (Price ~ Pimin) + fe(State)
result = reg(df, model, save = true)
@test residuals(result)[1:3] ≈ [ -16.925748, -14.835710, -12.017037] atol = 1e-4

#r2 with weights when saving residuals
m = @formula Sales ~ Price
result = reg(df, save = :residuals, weights = :Pop, m)
@test r2(result) ≈ 0.24654260 atol = 1e-4

# test different arguments for the keyword argument save
model = @formula Sales ~ Price + fe(State)
result = reg(df, model, save = true)
@test residuals(result) !== nothing
@test "fe_State" ∈ names(fe(result))

model = @formula Sales ~ Price + fe(State)
result = reg(df, model, save = :residuals)
@test residuals(result) !== nothing
@test "fe_State" ∉ names(fe(result))

model = @formula Sales ~ Price + fe(State)
result = reg(df, model, save = :fe)
@test residuals(result) === nothing
@test "fe_State" ∈ names(fe(result))

# iv recategorized
df.Pimin2 = df.Pimin
m = @formula Sales ~ (Pimin2 + Price ~ NDI + Pimin)
result = reg(df, m)
yhat = predict(result, df)
res = residuals(result, df)

m2 = @formula Sales ~ Pimin2 + (Price ~ NDI + Pimin)
result2 = reg(df, m2)
yhat2 = predict(result2, df)
res2 = residuals(result2, df)
@test yhat ≈ yhat2
@test res ≈ res2

m3 = @formula Sales ~ Pimin2 + (Price ~ NDI)
result3 = reg(df, m3)
yhat3 = predict(result3, df)
res3 = residuals(result3, df)
@test yhat ≈ yhat3
@test res ≈ res3

m4 = @formula Sales ~ (Price + Pimin2 ~ NDI + Pimin)
result4 = reg(df, m4)
yhat4 = predict(result4, df)
res4 = residuals(result4, df)
@test yhat ≈ yhat4
@test res ≈ res4

m5 = @formula Sales ~ (Price ~ NDI + Pimin) +  Pimin2
result5 = reg(df, m5)
yhat5 = predict(result5, df)
res5 = residuals(result5, df)
@test yhat ≈ yhat5
@test res ≈ res5

##############################################################################
##
## Saved FixedEffects
##
##############################################################################
# check save does not change r2
model1 = @formula Sales ~ Price
result1 = reg(df, model1, weights = :Pop)
model2 = @formula Sales ~ Price
result2 = reg(df, model2, weights = :Pop)
@test r2(result1) ≈ r2(result2)

model = @formula Sales ~ Price + fe(Year)
result = reg(df, model, save = true)
@test fe(result)[1, :fe_Year] ≈ 164.77833189721005
@test size(fe(result), 2) == 1
@test size(fe(result, keepkeys = true), 2) == 2

model = @formula Sales ~ Price + fe(Year) + fe(State)
result = reg(df, model, save = true)
@test fe(result)[1, :fe_Year] + fe(result)[1, :fe_State] ≈ 140.6852 atol = 1e-3

model = @formula Sales ~ Price + Year&fe(State)
result = reg(df, model, save = true)
@test fe(result)[1, Symbol("fe_State&Year")] ≈ 1.742779  atol = 1e-3

model = @formula Sales ~ Price + fe(State) + Year&fe(State)
result = reg(df, model, save = true)
@test fe(result)[1, :fe_State] ≈ -91.690635 atol = 1e-1

model = @formula Sales ~ Price + fe(State)
result = reg(df, model, subset = df.State .<= 30, save = true)
@test fe(result)[1, :fe_State] ≈  124.913976 atol = 1e-1
@test ismissing(fe(result)[1380 , :fe_State])

model = @formula Sales ~ Price + fe(Year)
result = reg(df, model, weights = :Pop, save = true)
@test fe(result)[2, :fe_Year] -  fe(result)[1, :fe_Year] ≈ -3.0347149502496222

# fixed effects
df.Price2 = df.Price
model = @formula Sales ~ Price + Price2 + fe(Year)
result = reg(df, model, save = true)
@test fe(result)[1, :fe_Year] ≈ 164.77833189721005

# iv
model = @formula Sales ~ (State ~ Price) + fe(Year)
result = reg(df, model, save = true)
@test fe(result)[1, :fe_Year] ≈ -167.48093490413623

# weights
model = @formula Sales ~ Price + fe(Year)
result = reg(df, model, weights = :Pop, save = true)
@test fe(result)[2, :fe_Year] -  fe(result)[1, :fe_Year] ≈ -3.0347149502496222

# IV and weights
model = @formula Sales ~ (Price ~ Pimin) + fe(Year)
result = reg(df, model, weights = :Pop, save = true)
@test fe(result)[1, :fe_Year] ≈ 168.24688 atol = 1e-4

# IV, weights and both year and state fixed effects
model = @formula Sales ~ (Price ~ Pimin) + fe(State) + fe(Year)
result = reg(df, model, weights = :Pop, save = true)
@test fe(result)[1, :fe_Year] + fe(result)[1, :fe_State]≈ 147.84145 atol = 1e-4

# subset with IV
model = @formula Sales ~ (Price ~ Pimin) + fe(Year)
result = reg(df, model, subset = df.State .<= 30, save = true)
@test fe(result)[1, :fe_Year] ≈ 164.05245824240276 atol = 1e-4
@test ismissing(fe(result)[811, :fe_Year])

# subset with IV, weights and year fixed effects
model = @formula Sales ~ (Price ~ Pimin) + fe(Year)
result = reg(df, model, subset = df.State .<= 30, weights = :Pop, save = true)
@test fe(result)[1, :fe_Year] ≈ 182.71915 atol = 1e-4

# subset with IV, weights and year fixed effects
model = @formula Sales ~ (Price ~ Pimin) + fe(State) + fe(Year)
result = reg(df, model, subset = df.State .<= 30, weights = :Pop, save = true)
@test fe(result)[1, :fe_Year] + fe(result)[1, :fe_State] ≈ 158.91798 atol = 1e-4

methods_vec = [:cpu]
if FixedEffectModels.FixedEffects.has_CUDA()
	push!(methods_vec, :gpu)
end
for method in methods_vec
	local model = @formula Sales ~ Price + fe(Year)
	local result = reg(df, model, save = true, method = method, double_precision = false)
	@test fe(result)[1, :fe_Year] ≈ 164.7 atol = 1e-1
end


text:
struct Row
    ind1::Int
    ind2::Int
    char::Char
    string::AbstractString
end

function torow(line::AbstractString)
    ind1, ind2, char, string = split(line, r"(:? )|-")
    ind1, ind2 = parse.(Int, [ind1, ind2])
    char = first(char)
    return Row(ind1, ind2, char, string)
end

function testline(line::AbstractString)
    row = torow(line)
    counts = count(i -> (i == row.char), row.string)
    return row.ind1 <= counts <= row.ind2
end

function testline2(line::AbstractString)
    row = torow(line)
    # Only one of the two must be equal to row.char
    return sum([row.string[row.ind1], row.string[row.ind2]] .== row.char) == 1
end

function part1()
    return sum(testline(line) for line in readlines("day02.txt"))
end

function part2()
    return sum(testline2(line) for line in readlines("day02.txt"))
end

println("Part1: $(part1())")
println("Part2: $(part2())")


text:
using Test

using MathOptInterface
const MOI = MathOptInterface
const MOIT = MathOptInterface.Test
const MOIU = MathOptInterface.Utilities
const MOIB = MathOptInterface.Bridges

include("../utilities.jl")

mock = MOIU.MockOptimizer(MOIU.UniversalFallback(MOIU.Model{Float64}()))
config = MOIT.TestConfig()

@testset "LogDet" begin
    bridged_mock = MOIB.Constraint.LogDet{Float64}(mock)

    @testset "logdet1test" begin
        var_primal = [0, 1, 0, 1, 1, 0, 1, 0, 0, 1]
        exp_duals = [[-1, -1, 1], [-1, -1, 1]]
        psd_dual = [1, 0, 1, -1, 0, 1, 0, -1, 0, 1]
        mock.optimize! =
            (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
                mock,
                var_primal,
                (MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) => [-1],
                (MOI.VectorAffineFunction{Float64}, MOI.Nonnegatives) =>
                    [[1, 1]],
                (MOI.VectorAffineFunction{Float64}, MOI.ExponentialCone) =>
                    exp_duals,
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                ) => [psd_dual],
            )

        MOIT.logdett1vtest(bridged_mock, config)
        MOIT.logdett1ftest(bridged_mock, config)

        # set primal/dual start is not yet implemented for LogDet bridge
        ci = first(
            MOI.get(
                bridged_mock,
                MOI.ListOfConstraintIndices{
                    MOI.VectorAffineFunction{Float64},
                    MOI.LogDetConeTriangle,
                }(),
            ),
        )
        test_delete_bridge(
            bridged_mock,
            ci,
            5,
            (
                (MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}, 0),
                (MOI.VectorAffineFunction{Float64}, MOI.ExponentialCone, 0),
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                    0,
                ),
            ),
        )
    end

    @testset "logdet2test" begin
        var_primal =
            [log(5), 3, 2, 2 / 3, 1, 1 / 3, 2.5, log(3), log(2 / 3), log(2.5)]
        exp_duals = [
            [-1, log(3) - 1, 1 / 3],
            [-1, log(2 / 3) - 1, 3 / 2],
            [-1, log(2.5) - 1, 0.4],
        ]
        psd_dual = [
            1,
            -1,
            1.6,
            0,
            -0.2,
            0.4,
            -1 / 3,
            0,
            0,
            1 / 3,
            1,
            -1.5,
            0,
            0,
            1.5,
            0,
            0.2,
            -0.4,
            0,
            0,
            0.4,
        ]
        mock.optimize! =
            (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
                mock,
                var_primal,
                (MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}) => [-1.0],
                (MOI.VectorAffineFunction{Float64}, MOI.ExponentialCone) =>
                    exp_duals,
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                ) => [psd_dual],
            )

        MOIT.logdett2test(bridged_mock, config)

        # set primal/dual start is not yet implemented for LogDet bridge
        ci = first(
            MOI.get(
                bridged_mock,
                MOI.ListOfConstraintIndices{
                    MOI.VectorAffineFunction{Float64},
                    MOI.LogDetConeTriangle,
                }(),
            ),
        )
        test_delete_bridge(
            bridged_mock,
            ci,
            1,
            (
                (MOI.ScalarAffineFunction{Float64}, MOI.LessThan{Float64}, 0),
                (MOI.VectorAffineFunction{Float64}, MOI.ExponentialCone, 0),
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                    0,
                ),
            ),
        )
    end
end

@testset "RootDet" begin
    bridged_mock = MOIB.Constraint.RootDet{Float64}(mock)

    @testset "rootdet1test" begin
        var_primal = [1, 1, 0, 1, 1, 0, 1]
        geomean_dual = [-1, 0.5, 0.5]
        psd_dual = [0.5, 0, 0.5, -0.5, 0, 0.5, 0, -0.5, 0, 0.5]
        mock.optimize! =
            (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
                mock,
                var_primal,
                (MOI.VectorAffineFunction{Float64}, MOI.Nonnegatives) =>
                    [[0.5, 0.5]],
                (MOI.VectorAffineFunction{Float64}, MOI.GeometricMeanCone) => [geomean_dual],
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                ) => [psd_dual],
            )

        MOIT.rootdett1vtest(bridged_mock, config)
        MOIT.rootdett1ftest(bridged_mock, config)

        # set primal/dual start is not yet implemented for RootDet bridge
        ci = first(
            MOI.get(
                bridged_mock,
                MOI.ListOfConstraintIndices{
                    MOI.VectorAffineFunction{Float64},
                    MOI.RootDetConeTriangle,
                }(),
            ),
        )
        test_delete_bridge(
            bridged_mock,
            ci,
            4,
            (
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.RotatedSecondOrderCone,
                    0,
                ),
                (MOI.VectorAffineFunction{Float64}, MOI.GeometricMeanCone, 0),
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                    0,
                ),
            ),
        )
    end

    @testset "rootdet2test" begin
        var_primal = [5^inv(3), 3, 2, 2 / 3, 1, 1 / 3, 2.5]
        multiplier = 5^inv(3) / 3
        geomean_dual = vcat(-1, multiplier * [inv(3), 1.5, 0.4])
        psd_dual =
            multiplier * [
                1,
                -1,
                1.6,
                0,
                -0.2,
                0.4,
                -1 / 3,
                0,
                0,
                1 / 3,
                1,
                -1.5,
                0,
                0,
                1.5,
                0,
                0.2,
                -0.4,
                0,
                0,
                0.4,
            ]
        mock.optimize! =
            (mock::MOIU.MockOptimizer) -> MOIU.mock_optimize!(
                mock,
                var_primal,
                (MOI.VectorAffineFunction{Float64}, MOI.GeometricMeanCone) => [geomean_dual],
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                ) => [psd_dual],
            )

        MOIT.rootdett2test(bridged_mock, config)

        # set primal/dual start is not yet implemented for RootDet bridge
        ci = first(
            MOI.get(
                bridged_mock,
                MOI.ListOfConstraintIndices{
                    MOI.VectorAffineFunction{Float64},
                    MOI.RootDetConeTriangle,
                }(),
            ),
        )
        test_delete_bridge(
            bridged_mock,
            ci,
            1,
            (
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.RotatedSecondOrderCone,
                    0,
                ),
                (MOI.VectorAffineFunction{Float64}, MOI.GeometricMeanCone, 0),
                (
                    MOI.VectorAffineFunction{Float64},
                    MOI.PositiveSemidefiniteConeTriangle,
                    0,
                ),
            ),
        )
    end
end


text:
# NOTE: This file will soon be automatically generated

import CLIMA

struct ParameterSet <: CLIMA.Parameters.AbstractParameterSet end

# Properties of dry air
CLIMA.Parameters.Planet.molmass_dryair(ps::ParameterSet) = 28.97e-3
CLIMA.Parameters.Planet.R_d(ps::ParameterSet) =
    gas_constant() / molmass_dryair(ps)
CLIMA.Parameters.Planet.kappa_d(ps::ParameterSet) = 2 / 7
CLIMA.Parameters.Planet.cp_d(ps::ParameterSet) = R_d(ps) / kappa_d(ps)
CLIMA.Parameters.Planet.cv_d(ps::ParameterSet) = cp_d(ps) - R_d(ps)

# Properties of water
CLIMA.Parameters.Planet.ρ_cloud_liq(ps::ParameterSet) = 1e3
CLIMA.Parameters.Planet.ρ_cloud_ice(ps::ParameterSet) = 916.7
CLIMA.Parameters.Planet.molmass_water(ps::ParameterSet) = 18.01528e-3
CLIMA.Parameters.Planet.molmass_ratio(ps::ParameterSet) =
    molmass_dryair(ps) / molmass_water(ps)
CLIMA.Parameters.Planet.R_v(ps::ParameterSet) =
    gas_constant() / molmass_water(ps)
CLIMA.Parameters.Planet.cp_v(ps::ParameterSet) = 1859
CLIMA.Parameters.Planet.cp_l(ps::ParameterSet) = 4181
CLIMA.Parameters.Planet.cp_i(ps::ParameterSet) = 2100
CLIMA.Parameters.Planet.cv_v(ps::ParameterSet) = cp_v(ps) - R_v(ps)
CLIMA.Parameters.Planet.cv_l(ps::ParameterSet) = cp_l(ps)
CLIMA.Parameters.Planet.cv_i(ps::ParameterSet) = cp_i(ps)
CLIMA.Parameters.Planet.T_freeze(ps::ParameterSet) = 273.15
CLIMA.Parameters.Planet.T_min(ps::ParameterSet) = 150.0
CLIMA.Parameters.Planet.T_max(ps::ParameterSet) = 1000.0
CLIMA.Parameters.Planet.T_icenuc(ps::ParameterSet) = 233.00
CLIMA.Parameters.Planet.T_triple(ps::ParameterSet) = 273.16
CLIMA.Parameters.Planet.T_0(ps::ParameterSet) = T_triple(ps)
CLIMA.Parameters.Planet.LH_v0(ps::ParameterSet) = 2.5008e6
CLIMA.Parameters.Planet.LH_s0(ps::ParameterSet) = 2.8344e6
CLIMA.Parameters.Planet.LH_f0(ps::ParameterSet) = LH_s0(ps) - LH_v0(ps)
CLIMA.Parameters.Planet.e_int_v0(ps::ParameterSet) =
    LH_v0(ps) - R_v(ps) * T_0(ps)
CLIMA.Parameters.Planet.e_int_i0(ps::ParameterSet) = LH_f0(ps)
CLIMA.Parameters.Planet.press_triple(ps::ParameterSet) = 611.657

# Properties of sea water
CLIMA.Parameters.Planet.ρ_ocean(ps::ParameterSet) = 1.035e3
CLIMA.Parameters.Planet.cp_ocean(ps::ParameterSet) = 3989.25

# Planetary parameters
CLIMA.Parameters.Planet.planet_radius(ps::ParameterSet) = 6.371e6
CLIMA.Parameters.Planet.day(ps::ParameterSet) = 86400
CLIMA.Parameters.Planet.Omega(ps::ParameterSet) = 7.2921159e-5
CLIMA.Parameters.Planet.grav(ps::ParameterSet) = 9.81
CLIMA.Parameters.Planet.year_anom(ps::ParameterSet) = 365.26 * day(ps)
CLIMA.Parameters.Planet.orbit_semimaj(ps::ParameterSet) = 1 * astro_unit()
CLIMA.Parameters.Planet.TSI(ps::ParameterSet) = 1362
CLIMA.Parameters.Planet.MSLP(ps::ParameterSet) = 1.01325e5

max_stars_repo_name: trontrytel/CLIMA

max_issues_repo_name: trontrytel/CLIMA

max_forks_repo_name: trontrytel/CLIMA


text:
using NoncommutativeGraphs, Subspaces
using Convex, SCS, LinearAlgebra
using Random, RandomMatrices
using Test

function random_bounded(n)
    U = rand(Haar(2), n)
    return Hermitian(U' * Diagonal(rand(n)) * U)
end

function basis_vec(dims::Tuple, idx::Tuple)
    m = zeros(dims)
    m[idx...] = 1
    return m
end

eye(n) = Matrix(1.0*I, (n,n))

solver_eps = 1e-7
tol = 1e-6

max_stars_repo_name: dstahlke/NoncommutativeGraphs

max_issues_repo_name: dstahlke/NoncommutativeGraphs

max_forks_repo_name: dstahlke/NoncommutativeGraphs


text:
import AstroTime.AstroDates: DateTime, date, time
import ERFA

@testset "Offsets" begin
    @testset "UTC" begin
        before_leap = UTCEpoch(2012, 6, 30, 23, 59, 59.0)
        before_leap_dt = DateTime(2012, 6, 30, 23, 59, 59.0)
        start_leap = UTCEpoch(2012, 6, 30, 23, 59, 60.0)
        start_leap_dt = DateTime(2012, 6, 30, 23, 59, 60.0)
        during_leap = UTCEpoch(2012, 6, 30, 23, 59, 60.3)
        during_leap_dt = DateTime(2012, 6, 30, 23, 59, 60.3)
        after_leap = UTCEpoch(2012, 7, 1)
        after_leap_dt = DateTime(2012, 7, 1)

        @test tai_offset(UTC, date(before_leap_dt), time(before_leap_dt)) == 34.0
        @test tai_offset(UTC, date(start_leap_dt), time(start_leap_dt)) == 34.0
        @test tai_offset(UTC, date(during_leap_dt), time(during_leap_dt)) == 34.0
        @test tai_offset(UTC, date(after_leap_dt), time(after_leap_dt)) == 35.0
    end
    @testset "TDB" begin
        ep = UTCEpoch(2000, 1, 1)
        @test TDBEpoch(ep) ≈ TDBEpoch(ep, 0.0, 0.0, 0.0) rtol=1e-3
        jd1, jd2 = value.(julian_twopart(ep))
        ut = fractionofday(UT1Epoch(ep))
        elong, u, v = abs.(randn(3)) * 1000
        exp = tai_offset(TT, ep) + ERFA.dtdb(jd1, jd2, ut, elong, u, v)
        act = tai_offset(TDB, ep, elong, u, v)
        @test act ≈ exp
    end
end


text:
#Inputs to the OperationModel

@testset "Operation Model Constructor with Params" begin
    #These tests with Unit Commitment must not add Parameters to the Canonical Model.
    # Once a formulation with Parameters is available the tests need to be updated
    #=
    devices = Dict{Symbol, PSI.DeviceModel}(:Generators => PSI.DeviceModel(PSY.ThermalGen, PSI.ThermalUnitCommitment),
                                            :Loads => PSI.DeviceModel(PSY.PowerLoad, PSI.StaticPowerLoad))
    branches = Dict{Symbol, PSI.DeviceModel}(:Lines => PSI.DeviceModel(PSY.Branch, PSI.SeriesLine))
    services = Dict{Symbol, PSI.ServiceModel}(:Reserves => PSI.ServiceModel(PSY.Reserve, PSI.AbstractReservesForm))
    PTDF, A = PowerSystems.buildptdf(branches5, nodes5)

    op_model = PSI.PowerOperationModel(TestOptModel, PM.StandardACPForm, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions) && :var_reactive in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PM.DCPlosslessForm, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.StandardPTDFForm, devices, branches, services, sys5b; PTDF = PTDF)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.CopperPlatePowerModel, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))
    =#

    devices = Dict{Symbol, PSI.DeviceModel}(:Generators => PSI.DeviceModel(PSY.ThermalGen, PSI.ThermalDispatch),
    :Loads => PSI.DeviceModel(PSY.PowerLoad, PSI.StaticPowerLoad))
    branches = Dict{Symbol, PSI.DeviceModel}(:Lines => PSI.DeviceModel(PSY.Branch, PSI.SeriesLine))
    services = Dict{Symbol, PSI.ServiceModel}(:Reserves => PSI.ServiceModel(PSY.Reserve, PSI.AbstractReservesForm))
    PTDF, A = PowerSystems.buildptdf(branches5, nodes5)

    op_model = PSI.PowerOperationModel(TestOptModel, PM.StandardACPForm, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions) && :var_reactive in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PM.DCPlosslessForm, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.StandardPTDFForm, devices, branches, services, sys5b; PTDF = PTDF)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.CopperPlatePowerModel, devices, branches, services, sys5b)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test (:params in keys(op_model.canonical_model.JuMPmodel.ext))
end

@testset "Operation Model Constructor without Params" begin
    devices = Dict{Symbol, PSI.DeviceModel}(:Generators => PSI.DeviceModel(PSY.ThermalGen, PSI.ThermalUnitCommitment),
    :Loads => PSI.DeviceModel(PSY.PowerLoad, PSI.StaticPowerLoad))
    branches = Dict{Symbol, PSI.DeviceModel}(:Lines => PSI.DeviceModel(PSY.Branch, PSI.SeriesLine))
    services = Dict{Symbol, PSI.ServiceModel}(:Reserves => PSI.ServiceModel(PSY.Reserve, PSI.AbstractReservesForm))
    PTDF, A = PowerSystems.buildptdf(branches5, nodes5)

    op_model = PSI.PowerOperationModel(TestOptModel, PM.StandardACPForm, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions) && :var_reactive in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PM.DCPlosslessForm, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.StandardPTDFForm, devices, branches, services, sys5b; PTDF = PTDF, parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.CopperPlatePowerModel, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    devices = Dict{Symbol, PSI.DeviceModel}(:Generators => PSI.DeviceModel(PSY.ThermalGen, PSI.ThermalDispatch),
                                            :Loads => PSI.DeviceModel(PSY.PowerLoad, PSI.StaticPowerLoad))
    branches = Dict{Symbol, PSI.DeviceModel}(:Lines => PSI.DeviceModel(PSY.Branch, PSI.SeriesLine))
    services = Dict{Symbol, PSI.ServiceModel}(:Reserves => PSI.ServiceModel(PSY.Reserve, PSI.AbstractReservesForm))
    PTDF, A = PowerSystems.buildptdf(branches5, nodes5)

    op_model = PSI.PowerOperationModel(TestOptModel, PM.StandardACPForm, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions) && :var_reactive in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PM.DCPlosslessForm, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.StandardPTDFForm, devices, branches, services, sys5b; PTDF = PTDF, parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

    op_model = PSI.PowerOperationModel(TestOptModel, PSI.CopperPlatePowerModel, devices, branches, services, sys5b; parameters = false)
    @test :var_active in keys(op_model.canonical_model.expressions)
    @test !(:params in keys(op_model.canonical_model.JuMPmodel.ext))

end

@testset "Build Operation Models" begin
    SCED = PSI.SCEconomicDispatch(sys5b; optimizer = GLPK_optimizer);
    OPF = PSI.OptimalPowerFlow(sys5b, PM.StandardACPForm, optimizer = ipopt_optimizer)
    UC = PSI.UnitCommitment(sys5b, PM.DCPlosslessForm; optimizer = GLPK_optimizer)

    ED_rts_p = PSI.EconomicDispatch(sys_rts, PSI.CopperPlatePowerModel; optimizer = GLPK_optimizer);
    ED_rts = PSI.EconomicDispatch(sys_rts, PSI.CopperPlatePowerModel; optimizer = GLPK_optimizer, parameters = false);
    # These other tests can be enabled when CDM parser get the correct HVDC type.
    #OPF_rts = PSI.OptimalPowerFlow(sys_rts, PSI.CopperPlatePowerModel, optimizer = ipopt_optimizer)
    #UC_rts = PSI.UnitCommitment(sys_rts, PSI.CopperPlatePowerModel; optimizer = GLPK_optimizer, parameters = false)
end


text:
#Code written by Skipper af Rosenborg and Esben Bager
#File with functions generating LASSO results
using JuMP
using Gurobi
env = Gurobi.Env()

function generateLassoModel(Xtrain, Ytrain, gamma)
	bCols = size(Xtrain)[2]
	M = JuMP.Model(solver = GurobiSolver(env, OutputFlag = 0, Threads=1))
	@variables M begin
			b[1:bCols]
			t[1:bCols]
			w
	end
	@objective(M,Min,w+gamma*ones(bCols)'*t)
	@constraint(M, soc, norm( [1-w; sqrt(2)*(Xtrain*b-Ytrain)] ) <= 1+w) #second order cone constraint
	@constraint(M,  b .<= t)
	@constraint(M, -t .<= b)

	solve(M)

#	println("solved by worker $(myid())")
	return getvalue(b)

#	return M
end

function lassoValidationAndTest(gamma, xTrain, yTrain, xVali, yVali, xTest, yTest)
	beta = generateLassoModel(xTrain, yTrain, gamma)

	#Shrink LASSO values
	for i in 1:length(beta)
		if beta[i] <= 1e-6
			if beta[i] >= -1e-6
				beta[i] = 0
			end
		end
	end

	k = countnz(beta)
	if k == 0
		maxPairwise = 0
	else
		maxPairwise = findMaxCor(xTrain, beta)
	end

	valiRsqr = getRSquared(xVali, yVali, beta)
	testRsqr = getRSquared(xTest, yTest, beta)

	return maxPairwise, testRsqr, valiRsqr, k, beta
end

function processOutput(Xtrain, Ytrain, Xpred, Ypred, bSolved)
	for i in 1:length(bSolved)
		if bSolved[i] <= 1e-6
			if bSolved[i] >= -1e-6
				bSolved[i] = 0
			end
		end
	end

	#In-Sample R-squared value
	errors = (Ytrain-Xtrain*bSolved)
	errorTotal = sum(errors[i]^2 for i=1:length(errors))
	errorsMean = Ytrain-mean(Ytrain)
	errorMeanTotal = sum(errorsMean[i]^2 for i=1:length(errorsMean))
	ISRsquared = 1-(errorTotal/errorMeanTotal)

	#Indicator Results
	YpredValue = Ypred[1]
	Yestimate = Xpred'*bSolved
	YestimateValue = Yestimate[1]
	if YestimateValue > 0
		Indicator = 1
	else
		Indicator = 0
	end

	return ISRsquared, Indicator, YestimateValue, bSolved
end

function generatSolveAndProcess(Xtrain, Ytrain, Xpred, Ypred, gamma)
	bSolved = generateLassoModel(Xtrain, Ytrain, gamma)

	ISRsquared, Indicator, YestimateValue, bSolved = processOutput(Xtrain, Ytrain, Xpred, Ypred, bSolved)

	return ISRsquared, Indicator, YestimateValue, bSolved
end

max_stars_repo_name: skipperafrosenborg/Thesis

max_issues_repo_path: ML/Code to hand in/ParallelModelGeneration.jl

max_issues_repo_name: skipperafrosenborg/Thesis

max_forks_repo_name: skipperafrosenborg/Thesis


text:
#=
These functions are taken from http://pordlabs.ucsd.edu/matlab/stream.htm
Kirill K. Pankratov, March 7, 1994.
=#

function cumsimp(vec::Vector{T}) where T
    c1 = 3//8
    c2 = 6//8
    c3 = -1//8

    # Interpolate values of Y to all midpoints
    f = fill(zero(T),length(vec))
    for i in 2:length(f)-1
        f[i] = c1*vec[i-1]+c2*vec[i]+c3*vec[i+1]
    end
    for i in 3:length(f)
        f[i] += c3*vec[i-2]+c2*vec[i-1]+c1*vec[i]
    end
    f[2] *= 2
    f[end] *= 2
    # Now Simpson (1,4,1) rule
    for i in 2:length(f)
        f[i] = 2f[i]+vec[i-1]+vec[i]
    end
    f[1:end] .= cumsum(f)/6  # Cumulative sum, 6 - denom. from the Simpson rule
end

function cumsimp(U::AbstractMatrix{T}) where T
    ret = similar(U)
    for i in 1:size(U)[2]
        ret[:,i] .= cumsimp(U[:,i])
    end
    ret
end

function flowfun(u,v)
    lx = size(u,2)  # Size of the velocity matrices
    ly = size(u,1)
    # Integrate velocity fields to get potential and streamfunction
    # Use Simpson rule summation
    cx = cumsimp(v[1,:])  # Compute x-integration constant
    cy = cumsimp(u[:,1])  # Compute y-integration constant
    ψ = -cumsimp(u)+repeat(cx',ly,1)
    ψ = (ψ+cumsimp(v')'-repeat(cy,1,lx))/2
end

# function cumsimp_naive(vec::Vector{T})::Vector{T} where T
#     c1 = 3//8
#     c2 = 6//8
#     c3 = -1//8
#
#     lv = length(vec)
#     num = 1:lv-2
#     vec_itp = fill(zero(T),length(vec))
#     # Interpolate values of Y to all midpoints
#     @. vec_itp[num+1] = c1*vec[num]+c2*vec[num+1]+c3*vec[num+2]
#     @. vec_itp[num+2] = vec_itp[num+2]+c3*vec[num]+c2*vec[num+1]+c1*vec[num+2]
#     vec_itp[2] *= 2
#     vec_itp[lv] *= 2
#     # Now Simpson (1,4,1) rule
#     @. vec_itp[num+1] = 2vec_itp[num+1]+vec[num]+vec[num+1]
#     vec_itp[:] .= cumsum(vec_itp)/6
# end