text:
"""
utility functions for Housedorff kernel
"""
module HFUtils
using ..CUDAGpuUtils, Logging,StaticArrays
 using CUDA, ..BasicStructs, Logging
 using ..CUDAGpuUtils ,..IterationUtils
export    clearLocArr,clearMainShmem,clearPadding,getIndexOfQueue
export @iter3dOuter,@iterDataBlock , calculateLoopsIter, @iterDataBlockZdeepest

"""
invoked before kernel execution in order to set number of needed loop iterations
dataBdim - dimensions of the data block
threadsXdim  - x dimension of the thread block
threadsYdim  - y dimension of the thread block
return tuple with numbers indicating how many iterations are needed in the loops of the kernel
"""
function calculateLoopsIter(dataBdim,threadsXdim,threadsYdim,metaDataDims,blocks)
    loopAXFixed= fld(dataBdim[2], threadsXdim)
    loopBXfixed= fld(dataBdim[3], threadsYdim)

    loopAYFixed= fld(dataBdim[1], threadsXdim)
    loopBYfixed= fld(dataBdim[3], threadsYdim)

    loopAZFixed= fld(dataBdim[1], threadsXdim)
    loopBZfixed= fld(dataBdim[2], threadsYdim)

    loopdataDimMainX = fld(dataBdim[1], threadsXdim)
    loopdataDimMainY = fld(dataBdim[2], threadsYdim)
    loopdataDimMainZ =dataBdim[2]
    inBlockLoopX,inBlockLoopY,inBlockLoopZ= (0,fld(dataBdim[2] ,threadsYdim),dataBdim[3]);
    metaDataLength= metaDataDims[1]*metaDataDims[2]*metaDataDims[3]

    loopMeta= fld(metaDataLength,blocks )
    loopWarpMeta= cld(metaDataLength,(blocks*threadsXdim ))
    resShmemTotalLength=(dataBdim[1]+2)*(dataBdim[2]+2)*(dataBdim[3]+2)
    sourceShmemTotalLength= dataBdim[1]*dataBdim[2]*dataBdim[3]
    clearIterResShmemLoop= fld(resShmemTotalLength,threadsXdim*threadsYdim)
    clearIterSourceShmemLoop= fld(sourceShmemTotalLength,threadsXdim*threadsYdim)

    shmemblockDataLenght = threadsXdim*threadsYdim*2
    shmemblockDataLoop = fld(shmemblockDataLenght,threadsXdim*threadsYdim)

    inBlockLoopXZIterWithPadding= fld(32,threadsYdim)

    return (inBlockLoopXZIterWithPadding,shmemblockDataLoop,shmemblockDataLenght,loopAXFixed,loopBXfixed,loopAYFixed,loopBYfixed,loopAZFixed,loopBZfixed,loopdataDimMainX,loopdataDimMainY,loopdataDimMainZ,inBlockLoopX,inBlockLoopY,inBlockLoopZ,metaDataLength,loopMeta,loopWarpMeta,clearIterResShmemLoop,clearIterSourceShmemLoop,resShmemTotalLength,sourceShmemTotalLength)
end

"""
1)   Left FP
2)   Left FN
3)   Right FP
4)   Right FN
5)   Posterior FP
6)   Posterior FN
7)   Anterior FP
8)   Anterior FN
9)   Top FP
10)   Top FN
11)   Bottom FP
12)   Bottom FN
13)   Total block Fp
14)   Total block Fn

xpos,ypos,zpos -current  position in x,y,z dimension
dataBdim - dimensions of the data block
boolSegm - true if we have true in algo arr
on the basis of the data it should give the index from 1 to 14 - to the appropriate queue
"""
function getIndexOfQueue(xpos,ypos,zpos, dataBdim,boolSegm)
    #we need to do so many != in order to deal with corners ...
    return (
     (xpos==1)*1
    +(xpos==dataBdim[1])*3
    +(ypos==1 && xpos!=1 && xpos!=dataBdim[1] )*5
    +(ypos==dataBdim[2] && xpos!=1 && xpos!=dataBdim[1]  )*7
    +(zpos==1  && xpos!=1 && xpos!=dataBdim[1]  && ypos!=1 && ypos!=dataBdim[2] )*9
    +(zpos==dataBdim[3] && xpos!=1 && xpos!=dataBdim[1]  && ypos!=1 && ypos!=dataBdim[2])*11
    +(xpos>1 && xpos<dataBdim[1] &&  ypos>1 && ypos<dataBdim[2] && zpos>1 && zpos<dataBdim[3])*13
    )+boolSegm# in that way we will get odd for fp an even for fn

end

"""
specialization of 2 dim iteration  for iterating over 3 dimensional metadata
    we join y and z iterations in order to increase occupancy for potentially small arrays
"""
macro iter3dOuter(metaDataDims,loopMeta,metaDataLength, ex)

    return  esc(quote
    linIdexMeta = UInt32(0)
    @unroll for j in 0:($loopMeta-1)
        linIdexMeta= blockIdxX()+ j*gridDimX()-1
        xMeta= rem(linIdexMeta,$metaDataDims[1])
        zMeta= fld((linIdexMeta),$metaDataDims[1]*$metaDataDims[2])
        yMeta= fld((linIdexMeta-((zMeta*$metaDataDims[1]*$metaDataDims[2] ) + xMeta )),$metaDataDims[1])
      $ex
    end
    linIdexMeta= blockIdxX()+ $loopMeta*gridDimX() -1
            if(linIdexMeta<$metaDataLength)
            xMeta= rem(linIdexMeta,$metaDataDims[1])
            zMeta= fld((linIdexMeta),$metaDataDims[1]*$metaDataDims[2])
            yMeta= fld((linIdexMeta-((zMeta*$metaDataDims[1]*$metaDataDims[2] ) + xMeta )),$metaDataDims[1])
        $ex
      end
    end)

end #iter3dOuter

"""
will enable iterating over the data of data block
"""
macro iterDataBlock(mainArrDims,dataBlockDims,loopXdim ,loopYdim,loopZdim,xMeta,yMeta,zMeta, ex)
    mainExp = generalizedItermultiDim(;arrDims=mainArrDims
    ,loopXdim
    ,loopYdim
    ,loopZdim
    # ,xCheck = :((xMeta* $dataBlockDims[1]+x)<=$mainArrDims[1] )
    # ,yCheck = :((yMeta* $dataBlockDims[2]+y)<=$mainArrDims[2])
    # ,zCheck = :( (zMeta* $dataBlockDims[3]+z)<=$mainArrDims[3])
    ,xCheck = :(((xdim * blockDimX())+threadIdxX()) <= $dataBlockDims[1] && x<=$mainArrDims[1])
    ,yCheck = :(((ydim * blockDimY())+threadIdxY()) <= $dataBlockDims[2] && y<=$mainArrDims[2])
    ,zCheck = :((zdim+1) <= $dataBlockDims[3]  &&   z <= $mainArrDims[3])
    ,zOffset= :($zMeta* ( ($dataBlockDims[3])))
    ,zAdd =:(zdim+1)
   ,yOffset = :(ydim* blockDimY()+$yMeta* $dataBlockDims[2])
   ,yAdd= :(threadIdxY())
   ,xOffset= :( (xdim * blockDimX()) +$xMeta* $dataBlockDims[1])
    ,xAdd= :(threadIdxX())
    ,isFullBoundaryCheckX =true
    , isFullBoundaryCheckY=true
    , isFullBoundaryCheckZ=true
    ,additionalActionBeforeX=  :(xpos= xdim * blockDimX()+threadIdxX() ;ypos= (ydim * blockDimY())+threadIdxY() ;zpos=(zdim+1))
    , ex = ex)
    return esc(:( $mainExp))
end

"""
will enable iterating over the data of data block
"""
macro iterDataBlockZdeepest(mainArrDims,dataBlockDims,loopXdim ,loopYdim,loopZdim,xMeta,yMeta,zMeta, ex,additionalActionAfterX)
    mainExp = generalizedItermultiDimZdeepest(;arrDims=mainArrDims
    ,loopXdim
    ,loopYdim
    ,loopZdim
    # ,xCheck = :((xMeta* $dataBlockDims[1]+x)<=$mainArrDims[1] )
    # ,yCheck = :((yMeta* $dataBlockDims[2]+y)<=$mainArrDims[2])
    # ,zCheck = :( (zMeta* $dataBlockDims[3]+z)<=$mainArrDims[3])
    ,xCheck = :(((xdim * blockDimX())+threadIdxX()) <= $dataBlockDims[1] && x<=$mainArrDims[1])
    ,yCheck = :(((ydim * blockDimY())+threadIdxY()) <= $dataBlockDims[2] && y<=$mainArrDims[2])
    ,zCheck = :((zdim+1) <= $dataBlockDims[3]  &&   z <= $mainArrDims[3])
    ,zOffset= :($zMeta* ( ($dataBlockDims[3])))
    ,zAdd =:(zdim+1)
   ,yOffset = :(ydim* blockDimY()+$yMeta* $dataBlockDims[2])
   ,yAdd= :(threadIdxY())
   ,xOffset= :( (xdim * blockDimX()) +$xMeta* $dataBlockDims[1])
    ,xAdd= :(threadIdxX())
    ,isFullBoundaryCheckX =true
    , isFullBoundaryCheckY=true
    , isFullBoundaryCheckZ=true
    ,additionalActionBeforeX=  :(xpos= xdim * blockDimX()+threadIdxX() ;ypos= (ydim * blockDimY())+threadIdxY())
    ,additionalActionBeforeZ=  :(zpos=(zdim+1))
    ,additionalActionAfterX=additionalActionAfterX
    , ex = ex)
    return esc(:( $mainExp))
end

#macro iterDataBlock(mainArrDims,dataBlockDims,loopXdim ,loopYdim,loopZdim,ex)
    #     mainExp = generalizedItermultiDim(;arrDims=mainArrDims
    #     ,loopXdim
    #     ,loopYdim
    #     ,loopZdim
    #     # ,xCheck = :((xMeta* $dataBlockDims[1]+x)<=$mainArrDims[1] )
    #     # ,yCheck = :((yMeta* $dataBlockDims[2]+y)<=$mainArrDims[2])
    #     # ,zCheck = :( (zMeta* $dataBlockDims[3]+z)<=$mainArrDims[3])
    #     ,xCheck = :(((xdim * blockDimX())+threadIdxX()) <= $dataBlockDims[1] && x<=$mainArrDims[1])
    #     ,yCheck = :(((ydim * blockDimY())+threadIdxY()) <= $dataBlockDims[2] && y<=$mainArrDims[2])
    #     ,zCheck = :((zdim+1) <= $dataBlockDims[3]  &&   z <= $mainArrDims[3])
    #     ,zOffset= :((zMeta-1)* ( ($dataBlockDims[3])))
    #     ,zAdd =:(zdim+1)
    #    ,yOffset = :(ydim* blockDimY()+(yMeta-1)* $dataBlockDims[2])
    #    ,yAdd= :(threadIdxY())
    #    ,xOffset= :( (xdim * blockDimX()) +(xMeta-1)* $dataBlockDims[1])
    #     ,xAdd= :(threadIdxX())
    #     ,isFullBoundaryCheckX =true
    #     , isFullBoundaryCheckY=true
    #     , isFullBoundaryCheckZ=true
    #     ,additionalActionBeforeX=  :(xpos= xdim * blockDimX()+threadIdxX() ;ypos= (ydim * blockDimY())+threadIdxY() ;zpos=(zdim+1))
    #     , ex = ex)
    #     return esc(:( $mainExp))
    # end

"""
this macro will be used to
    - decide wheather the part of the code will be executed at all - in order to avoid branching at runtime
    - select at compile time appropriate target for the statement (generally the global memory array )
    - all will be based on the  ConfigurtationStruct struct - conf
    entryName - will tell us what particularly we should look for in the conf
    if the entry of name entryName in conf will be true  we will return expression otherwise we will return nothing
    we will
    metricsArr - will be used from caller scope - and it is an array in global memory that will be
"""
macro isInConf(conf, entryName, ex)

    index = 0
    if(conf.sliceWiseMatrics)
        index=3# 3 becouse one will be added in a moment
    end
    # we look for the all of the properties that are true - we need to also remember that in case of sliceWiseMatrics - it will mark that we need 4 arrays
    for i in propertynames(eval(conf))
        if getproperty(conf, i)
            index+=1
        end#if
        if i == Symbol(entryName)
            break
        end
    end#for

    if getproperty(conf, Symbol(entryName))

    return esc(quote

        @inbounds metricsArr[$index]  = $ex # @inbounds fp[]+=

    end#quote
    )
end#if

end

end#HFUtils

# """
# we use it to access pading in the shmem
# change the indexing so  instead of transverse plane of threads we will get coronal
#     so we supply idx and idx y earlier z was constant and x,y changed so we considered top with z =1 and bottom z = 34
#     now we think x and y as would be displayed on transverse plane so in order to get anterior and posterior padding plane  we need to keep y as contant 1 or 33 and rest changing
#     and in order to get left and right we need to  keep x as contant to 1 or 34 (34 in case we have size of data block = 32)
#     so we will access shmem in changed indexing
# """
# function transvarseToCoronalThreadPlane(constantNumb::UInt8,shmem)::Bool
#     return shmem[threadIdxX(),constantNumb, threadIdxY()]
# end#transvarseToCoronalThreadPlane

# """
# we use it to access pading in the shmem
# change the indexing so  instead of transverse plane of threads we will get coronal
#     so we supply idx and idx y earlier z was constant and x,y changed so we considered top with z =1 and bottom z = 34
#     now we think x and y as would be displayed on transverse plane so in order to get anterior and posterior padding plane  we need to keep y as contant 1 or 33 and rest changing
#     and in order to get left and right we need to  keep x as contant to 1 or 34 (34 in case we have size of data block = 32)
#     so we will access shmem in changed indexing
# """
# function transvarseToSaggitalhreadPlane(constantNumb::UInt8,shmem)::Bool
#     return shmem[constantNumb,threadIdxX(), threadIdxY()]
# end#transvarseToSaggitalhreadPlane

# """
# clear main part of shared memory - padding will be cleared separately
# """
# function clearMainShmem(shmem)
#     @unroll for zIter in UInt8(1):UInt8(32)# most outer loop is responsible for z dimension
#         shmem[threadIdxX()+1,threadIdxY()+1,zIter+1]=0
#     end#for
# end#clearMainShmem

# """
# clear source shmem in shared memory """
# function clearSourceShmem(shmem)
#     @unroll for zIter in UInt8(1):UInt8(32)
#          shmem[threadIdxX(),threadIdxY(),zIter]=0
#     end#for
# end#clearMainShmem

# """
# set padding planes to 0
# """
# function clearPadding(shmem)
#     clearHalfOfPadding(shmem,UInt8(1))
#     clearHalfOfPadding(shmem,UInt8(34))
# end#clearPadding
# """
# helper for clearPadding
# """
# function clearHalfOfPadding(shmem,constantNumb::UInt8)
#     shmem[constantNumb,threadIdxX()+1, threadIdxY()+1]=false
#     shmem[threadIdxX()+1,constantNumb, threadIdxY()+1]=false
#     shmem[threadIdxX()+1,threadIdxY()+1, constantNumb]=false
# end


text:
module Constraints

export AbstractConstraint, FixedDistanceConstraint

abstract type AbstractConstraint end

struct FixedDistanceConstraint <: AbstractConstraint
    i::Int
    j::Int
    distance::Float64
end

end

ext: jl


text:
using Transpiler
using Base.Test

# TODO, GPUArrays tests most of the functionality of Transpiler.
# Maybe we should just depend on GPUArrays and run the OpenCL tests here?!

function test{T}(a::T, b)
    x = sqrt(sin(a) * b) / T(10.0)
    y = T(33.0)x + cos(b)
    y * T(10.0)
end

function fortest(x)
    acc = x
    for i = 1:5
        if i == 1
            acc += x
        elseif i == 2
            acc -= x
        else
            acc += x * x
        end
    end
    return acc
end

@testset "OpenCL Transpiler" begin
    include("opencl_funcs.jl")
end
# @testset "OpenGL Transpiler" begin
#     include("opengl_funcs.jl")
# end


text:
module ComposeBenchs

using InplaceArrays.Arrays
using InplaceArrays.Fields
import InplaceArrays.Fields: ∇
using InplaceArrays.Fields: MockField, MockBasis
using InplaceArrays.TensorValues
using FillArrays

@inline function repeat(n,f,args...)
  for i in 1:n
    f(args...)
  end
  nothing
end

@inline function loop(a,cache)
  for i in eachindex(a)
    ai = getindex!(cache,a,i)
  end
end

@noinline function loop_and_evaluate(ca,cai,cx,a,x)
  for i in eachindex(a)
    ai = getindex!(ca,a,i)
    xi = getindex!(cx,x,i)
    vi = evaluate!(cai,ai,xi)
  end
end

fun2(x,y) = 2*x
∇fun2(x,y) = VectorValue(2*one(x[1]),2*one(x[1]))
∇(::typeof(fun2)) = ∇fun2

function bench1(n)
  v = 3.0
  d = 2
  ndof = 8
  f = MockField{d}(v)
  b = MockBasis{d}(v,ndof)
  g = compose(fun2,f,v)
  np = 4
  p = Point(1,2)
  x = fill(p,np)
  cg = field_cache(g,x)
  @time repeat(n,evaluate!,cg,g,x)
  h = compose(fun2,f,f)
  ch = field_cache(h,x)
  @time repeat(n,evaluate!,ch,h,x)
  j = compose(fun2,b,f)
  cj = field_cache(j,x)
  @time repeat(n,evaluate!,cj,j,x)
end

function bench2(n)
  v = 3.0
  d = 2
  ndof = 8
  f = MockField{d}(v)
  b = MockBasis{d}(v,ndof)
  _g = compose(fun2,f,v)
  g = gradient(_g)
  np = 4
  p = Point(1,2)
  x = fill(p,np)
  cg = field_cache(g,x)
  @time repeat(n,evaluate!,cg,g,x)
  _h = compose(fun2,f,f)
  h = gradient(_h)
  ch = field_cache(h,x)
  @time repeat(n,evaluate!,ch,h,x)
  _j = compose(fun2,b,f)
  j = gradient(_j)
  cj = field_cache(j,x)
  @time repeat(n,evaluate!,cj,j,x)
end

function bench3(n)
  v = 3.0
  d = 2
  ndof = 8
  np = 4
  p = Point(1,2)
  x = fill(p,np)
  f = MockField{d}(v)
  b = MockBasis{d}(v,ndof)
  af = fill(f,n)
  ab = Fill(f,n)
  av = Fill(v,n)
  ax = fill(x,n)
  ag = compose(fun2,af,av)
  cag, cgi, cax = field_array_cache(ag,ax)
  @time loop_and_evaluate(cag,cgi,cax,ag,ax)
  ah = compose(fun2,af,ab)
  cah, chi, _ = field_array_cache(ah,ax)
  @time loop_and_evaluate(cah,chi,cax,ah,ax)
  aj = compose(fun2,av,ab)
  caj, cji, _ = field_array_cache(aj,ax)
  @time loop_and_evaluate(caj,cji,cax,aj,ax)
end

function bench4(n)
  v = 3.0
  d = 2
  ndof = 8
  np = 4
  p = Point(1,2)
  x = fill(p,np)
  f = MockField{d}(v)
  b = MockBasis{d}(v,ndof)
  af = fill(f,n)
  ab = Fill(f,n)
  av = Fill(v,n)
  ax = fill(x,n)
  _ag = compose(fun2,af,av)
  ag = gradient(_ag)
  cag, cgi, cax = field_array_cache(ag,ax)
  @time loop_and_evaluate(cag,cgi,cax,ag,ax)
  _ah = compose(fun2,af,ab)
  ah = gradient(_ah)
  cah, chi, _ = field_array_cache(ah,ax)
  @time loop_and_evaluate(cah,chi,cax,ah,ax)
  _aj = compose(fun2,av,ab)
  aj = gradient(_aj)
  caj, cji, _ = field_array_cache(aj,ax)
  @time loop_and_evaluate(caj,cji,cax,aj,ax)
end

foo(x) = 2*x[1]
∇foo(x) = VectorValue(2.0,0.0)
∇(::typeof(foo)) = ∇foo

function bench5(n)
  np = 4
  p = Point(1,2)
  x = fill(p,np)
  v = 3.0
  d = 2
  ndof = 8
  wi = 3.0
  w = fill(wi,ndof)
  l = n
  f = MockBasis{d}(v,ndof)
  af = Fill(f,l)
  ax = fill(x,l)
  aw = fill(w,l)
  _ag = lincomb(af,aw)
  ag = compose(foo,_ag)
  cag = array_cache(ag)
  #TODO
  #@time loop(ag,cag)
  agx = evaluate(ag,ax)
  cagx = array_cache(agx)
  @time loop(agx,cagx)
  g = testitem(ag)
  cg = field_cache(g,x)
  cax = array_cache(ax)
  #TODO
  #@time loop_and_evaluate(cag,cg,cax,ag,ax)
end

for n in (1,1,10,1000,100000)
  @eval begin
    println("+++ runing suite for n = $($n) +++")
    bench1($n)
    bench2($n)
    bench3($n)
    bench4($n)
    bench5($n)
  end
end

end # module


text:
# Device Control API

function deviceOpen(idxDevice::Int32)
  pHdwf = Ref{Hdwf}(0)
  succ = ccall((:FDwfDeviceOpen,libdwf),Cint,(Cint,Ptr{Hdwf}),
    idxDevice,pHdwf)
  succ == 0 && error("Error calling FDwfDeviceOpen.")
  return pHdwf[]
end

function deviceConfigOpen(idxDevice::Int32,idxCfg::Int32)
  pHdwf = Ref{Hdwf}(0)
  succ = ccall((:FDwfDeviceConfigOpen,libdwf),Cint,(Cint,Cint,Ptr{Hdwf}),
    idxDevice,idxCfg,pHdwf)
  succ == 0 && error("Error calling FDwfDeviceConfigOpen.")
  return pHdwf[]
end

function deviceClose(hdwf::Hdwf)
  succ = ccall((:FDwfDeviceClose,libdwf),Cint,(Hdwf,),
    hdwf)
  succ == 0 && error("Error calling FDwfDeviceClose.")
  return nothing
end

function deviceCloseAll()
  succ = ccall((:FDwfDeviceCloseAll,libdwf),Cint,()
    )
  succ == 0 && error("Error calling FDwfDeviceCloseAll.")
  return nothing
end

#### TODO ####
# DWFAPI BOOL FDwfDeviceAutoConfigureSet(HDWF hdwf, BOOL fAutoConfigure);
# DWFAPI BOOL FDwfDeviceAutoConfigureGet(HDWF hdwf, BOOL *pfAutoConfigure);
# DWFAPI BOOL FDwfDeviceReset(HDWF hdwf);
# DWFAPI BOOL FDwfDeviceEnableSet(HDWF hdwf, BOOL fEnable);
# DWFAPI BOOL FDwfDeviceTriggerInfo(HDWF hdwf, int *pfstrigsrc); // use IsBitSet
# DWFAPI BOOL FDwfDeviceTriggerSet(HDWF hdwf, int idxPin, TRIGSRC trigsrc);
# DWFAPI BOOL FDwfDeviceTriggerGet(HDWF hdwf, int idxPin, TRIGSRC *ptrigsrc);
# DWFAPI BOOL FDwfDeviceTriggerPC(HDWF hdwf);


text:
"""
    parserl(T, io)

Reads a line from `io` and parses the result to type `T`. Typically used for single-value lines.
"""
parserl(T, io) = parse(T, readline(io))

"""
    parsesrl(T, io; sep = ",", headskip = 0, tailskip = 0)

Reads a line from `io`, splits the line over `sep`. Use keywords `headskip` and `tailskip` to optionally exclude the head and tail of the string before splitting.

#Usage
julia> parsesrl(Int, io) # read "1,0,1"

"""
function parsesrl(T, io; sep = ",", headskip = 0, tailskip = 0)
    l = readline(io)
    i = 1 + headskip
    j = length(l) - tailskip
    l = l[i:j]
    parse.(T, split(l, sep))
end


text:
typealias Doc Vector{Vector{Token}}

function train!(::Type{Perceptron}, parser::DepParser, trainsents::Doc,
    testsents::Doc=Vector{Token}[]; beamsize=10, iter=20, progbar=true)
    # initialize parser.model
    initmodel!(parser, Perceptron)

    info("LOADING SENTENCES")
    info("WILL RUN $iter ITERATIONS")
    for i = 1:iter
        info("ITER $i TRAINING")
        progbar && ( p = Progress(length(trainsents), 1, "", 50) )
        res = map(trainsents) do s
            progbar && next!(p)
            s = State(s, parser)
            gold = beamsearch(s, 1, expandgold)
            pred = beamsearch(s, beamsize, expandpred)
            update!(parser.model, gold, pred)
            pred[end][1]
        end
        evaluate(parser, res)

        if !isempty(testsents)
            info("ITER $i TESTING")
            res = decode(Perceptron, parser, testsents)
            evaluate(parser, res)
        end
    end
end

function decode(::Type{Perceptron}, parser::DepParser, sents::Doc;
    beamsize=10, progbar=true)
    progbar && ( p = Progress(length(sents), 1, "", 50) )
    map(sents) do s
        progbar && next!(p)
        pred = State(s, parser)
        beamsearch(pred, beamsize, expandpred)[end][1]
    end
end

max_stars_repo_name: hshindo/JukaiNLP

max_issues_repo_name: hshindo/JukaiNLP

max_forks_repo_name: hshindo/JukaiNLP


text:
"""
Module này thực hiện tính toán với hàng đợi

# Cách dùng

- `m = TênMôHình(tham_số...)`: Tạo một mô hình, danh sách mô hình xem ở dưới
- `λe(m)`: Tính cường độ đến hiệu quả ``\\lambda_e``.
- `pn(m)`: Tính xác suất ``P\\left[N(t) = n\\right]``.
- `L(m)`: Lượng khách trung bình trong hệ thống
- `Lq(m)`: Lượng khách trung bình trong hàng đợi
- `W(m)`: Thời gian trung bình trong hệ thống
- `Wq(m)`: Thời gian trung bình trong hàng đợi

## Ví dụ
~~~julia
julia> include("$(@__FILE__)")
Main.MMCKModels

julia> using .MMCKModels

julia> m = MM1(6, 8)
MM1(6, 8, 0.75)

julia> L(m)
3.0

julia> Lq(m)
2.25

julia> W(m)
0.5
~~~

Các loại mô hình hàng đợi:
- Hiện có: `MM1`, `MM1K`, `MMC`, `MMCK`

Viết `?TênMôHình` để xem cách dùng.
"""
module MMs

#export MM1, MM1K, MMC,
#
export AbstractMMCK, MM1, MM1K, MMC, MMCK,
	L, Lq, W, Wq, λe, pn, cost

abstract type AbstractMMCK end;

include(joinpath(@__DIR__, "mm1.jl"))
include(joinpath(@__DIR__, "mm1k.jl"))
include(joinpath(@__DIR__, "mmc.jl"))
include(joinpath(@__DIR__, "mmck.jl"))

function λe(m::AbstractMMCK)
    hasproperty(m, :k) ? m.λ * (1 - pn(m, m.k)) : m.λ
end

end # module


text:
### A Pluto.jl notebook ###
# v0.17.5

using Markdown
using InteractiveUtils

# This Pluto notebook uses @bind for interactivity. When running this notebook outside of Pluto, the following 'mock version' of @bind gives bound variables a default value (instead of an error).
macro bind(def, element)
    quote
        local iv = try Base.loaded_modules[Base.PkgId(Base.UUID("6e696c72-6542-2067-7265-42206c756150"), "AbstractPlutoDingetjes")].Bonds.initial_value catch; b -> missing; end
        local el = $(esc(element))
        global $(esc(def)) = Core.applicable(Base.get, el) ? Base.get(el) : iv(el)
        el
    end
end

# ╔═╡ b65c67ec-b79f-4f0e-85e6-78ff22b279d4
using HypertextLiteral

# ╔═╡ a8c1e0d2-3604-4e1d-a87c-c8f5b86b79ed
md"""
# MultiCheckBox
"""

# ╔═╡ 1c914889-9c98-4aa5-b24f-3027b74feb4a
animals_count = Ref(0)

# ╔═╡ c8350f43-0d30-45d0-873b-ff56c5801ac1
md"""
## Definition
"""

# ╔═╡ 79c4dc76-efa2-4b7f-99ea-243e3a17f81c
function skip_as_script(m::Module)
	if isdefined(m, :PlutoForceDisplay)
		return m.PlutoForceDisplay
	else
		isdefined(m, :PlutoRunner) && parentmodule(m) == Main
	end
end

# ╔═╡ 499ca710-1a50-4aa1-87d8-d213416e8e30
if skip_as_script(@__MODULE__)
	import Pkg
	Pkg.activate(Base.current_project(@__DIR__))
	Text("Project env active")
end

# ╔═╡ 631c14bf-e2d3-4a24-8ddc-095a3dab80ef
import AbstractPlutoDingetjes.Bonds

# ╔═╡ c38de38d-e900-4309-a9f6-1392af2f245b
subarrays(x) = (
	x[collect(I)]
	for I in Iterators.product(Iterators.repeated([true,false],length(x))...) |> collect |> vec
)

# ╔═╡ 430e2c1a-832f-11eb-024a-13e3989fd7c2
begin
	export MultiCheckBox
    local result = begin
    """
    ```julia
    MultiCheckBox(options::Vector; [default::Vector], [orientation ∈ [:row, :column]], [select_all::Bool])
    ```

    A group of checkboxes - the user can choose which of the `options` to return.
    The value returned via `@bind` is a list containing the currently checked items.

    See also: [`MultiSelect`](@ref).

    `options` can also be an array of pairs `key::Any => value::String`. The `key` is returned via `@bind`; the `value` is shown.

    # Keyword arguments
    - `defaults` specifies which options should be checked initally.
    - `orientation` specifies whether the options should be arranged in `:row`'s `:column`'s.
    - `select_all` specifies whether or not to include a "Select All" checkbox.

    # Examples
    ```julia
    @bind snacks MultiCheckBox(["🥕", "🐟", "🍌"]))

    if "🥕" ∈ snacks
        "Yum yum!"
    end
    ```

    ```julia
    @bind functions MultiCheckBox([sin, cos, tan])

    [f(0.5) for f in functions]
    ```

    ```julia
    @bind snacks MultiCheckBox(["🥕" => "🐰", "🐟" => "🐱", "🍌" => "🐵"]; default=["🥕", "🍌"])
    ```

    ```julia
    @bind animals MultiCheckBox(["🐰", "🐱" , "🐵", "🐘", "🦝", "🐿" , "🐝",  "🐪"]; orientation=:column, select_all=true)
    ```
    """
    struct MultiCheckBox{BT,DT}
        options::AbstractVector{Pair{BT,DT}}
        default::Union{Missing,AbstractVector{BT}}
        orientation::Symbol
        select_all::Bool
    end
    end

    MultiCheckBox(options::AbstractVector{<:Pair{BT,DT}}; default=missing, orientation=:row, select_all=false) where {BT,DT} = MultiCheckBox(options, default, orientation, select_all)

    MultiCheckBox(options::AbstractVector{BT}; default=missing, orientation=:row, select_all=false) where BT = MultiCheckBox{BT,BT}(Pair{BT,BT}[o => o for o in options], default, orientation, select_all)

    function Base.show(io::IO, m::MIME"text/html", mc::MultiCheckBox)
        @assert mc.orientation == :column || mc.orientation == :row "Invalid orientation $(mc.orientation). Orientation should be :row or :column"

        defaults = coalesce(mc.default, [])

		# Old:
		# checked = [k in defaults for (k,v) in mc.options]
		#
		# More complicated to fix https://github.com/JuliaPluto/PlutoUI.jl/issues/106
		defaults_copy = copy(defaults)
		checked = [
			let
				i = findfirst(isequal(k), defaults_copy)
				if i === nothing
					false
				else
					deleteat!(defaults_copy, i)
					true
				end
			end
		for (k,v) in mc.options]

        show(io, m, @htl("""
        <plj-multi-checkbox style="flex-direction: $(mc.orientation);"></plj-multi-checkbox>
        <script type="text/javascript">
		const labels = $([string(v) for (k,v) in mc.options]);
		const values = $(1:length(mc.options));
		const checked = $(checked);
		const includeSelectAll = $(mc.select_all);

		const container = (currentScript ? currentScript : this.currentScript).previousElementSibling

		const my_id = crypto.getRandomValues(new Uint32Array(1))[0].toString(36)

		// Add checkboxes
		const inputEls = []
		for (let i = 0; i < labels.length; i++) {
			const boxId = `\${my_id}-box-\${i}`

			const item = document.createElement('div')

			const checkbox = document.createElement('input')
			checkbox.type = 'checkbox'
			checkbox.id = boxId
			checkbox.name = labels[i]
			checkbox.value = values[i]
			checkbox.checked = checked[i]
			inputEls.push(checkbox)
			item.appendChild(checkbox)

			const label = document.createElement('label')
			label.htmlFor = boxId
			label.innerText = labels[i]
			item.appendChild(label)

			container.appendChild(item)
		}

		function setValue() {
			container.value = inputEls.filter((o) => o.checked).map((o) => o.value)
		}
		// Add listeners
		function sendEvent() {
			setValue()
			container.dispatchEvent(new CustomEvent('input'))
		}

		function updateSelectAll() {}

		if (includeSelectAll) {
			// Add select-all checkbox.
			const selectAllItem = document.createElement('div')
			selectAllItem.classList.add(`select-all`)

			const selectID = `\${my_id}-select-all`

			const selectAllInput = document.createElement('input')
			selectAllInput.type = 'checkbox'
			selectAllInput.id = selectID
			selectAllItem.appendChild(selectAllInput)

			const selectAllLabel = document.createElement('label')
			selectAllLabel.htmlFor = selectID
			selectAllLabel.innerText = 'Select All'
			selectAllItem.appendChild(selectAllLabel)

			container.prepend(selectAllItem)

			function onSelectAllClick(event) {
				event.stopPropagation()
				inputEls.forEach((o) => (o.checked = this.checked))
				sendEvent()
			}
			selectAllInput.addEventListener('click', onSelectAllClick)
            selectAllInput.addEventListener('input', e => e.stopPropagation())

			/// Taken from: https://stackoverflow.com/questions/10099158/how-to-deal-with-browser-differences-with-indeterminate-checkbox
			/// Determine the checked state to give to a checkbox
			/// with indeterminate state, so that it becomes checked
			/// on click on IE, Chrome and Firefox 5+
			function getCheckedStateForIndeterminate() {
				// Create a unchecked checkbox with indeterminate state
				const test = document.createElement('input')
				test.type = 'checkbox'
				test.checked = false
				test.indeterminate = true

				// Try to click the checkbox
				const body = document.body
				body.appendChild(test) // Required to work on FF
				test.click()
				body.removeChild(test) // Required to work on FF

				// Check if the checkbox is now checked and cache the result
				if (test.checked) {
					getCheckedStateForIndeterminate = function () {
						return false
					}
					return false
				} else {
					getCheckedStateForIndeterminate = function () {
						return true
					}
					return true
				}
			}

			updateSelectAll = function () {
				const checked = inputEls.map((o) => o.checked)
				if (checked.every((x) => x)) {
					selectAllInput.checked = true
					selectAllInput.indeterminate = false
				} else if (checked.some((x) => x)) {
					selectAllInput.checked = getCheckedStateForIndeterminate()
					selectAllInput.indeterminate = true
				} else {
					selectAllInput.checked = false
					selectAllInput.indeterminate = false
				}
			}
			// Call once at the beginning to initialize.
			updateSelectAll()
		}

		function onItemClick(event) {
			event.stopPropagation()
			updateSelectAll()
			sendEvent()
		}
		setValue()
		inputEls.forEach((el) => el.addEventListener('click', onItemClick))
		inputEls.forEach((el) => el.addEventListener('input', e => e.stopPropagation()))

        </script>
        <style type="text/css">
		plj-multi-checkbox {
			display: flex;
			flex-wrap: wrap;
			/* max-height: 8em; */
		}

		plj-multi-checkbox * {
			display: flex;
		}

		plj-multi-checkbox > div {
			margin: 0.1em 0.3em;
			align-items: center;
		}

		plj-multi-checkbox label,
		plj-multi-checkbox input {
			cursor: pointer;
		}

		plj-multi-checkbox .select-all {
			font-style: italic;
			color: hsl(0, 0%, 25%, 0.7);
		}
		</style>
        """))
    end

    Base.get(select::MultiCheckBox) = Bonds.initial_value(select)
    Bonds.initial_value(select::MultiCheckBox{BT,DT}) where {BT,DT} =
        ismissing(select.default) ? BT[] : select.default
    Bonds.possible_values(select::MultiCheckBox) =
        subarrays((string(i) for i in 1:length(select.options)))

    function Bonds.transform_value(select::MultiCheckBox{BT,DT}, val_from_js) where {BT,DT}
        # val_from_js will be a vector of Strings, but let's allow Integers as well, there's no harm in that
        @assert val_from_js isa Vector

        val_nums = (
            v isa Integer ? v : tryparse(Int64, v)
            for v in val_from_js
        )

        BT[select.options[v].first for v in val_nums]
    end

    function Bonds.validate_value(select::MultiCheckBox, val)
        val isa Vector && all(val_from_js) do v
            val_num = v isa Integer ? v : tryparse(Int64, v)
            1 ≤ val_num ≤ length(select.options)
        end
    end
    result
end

# ╔═╡ 8bfaf4c8-557d-433e-a228-aac493746efc
@bind animals MultiCheckBox(["🐰", "🐱" , "🐵", "🐘", "🦝", "🐿" , "🐝",  "🐪"]; orientation=:column, select_all=true)

# ╔═╡ 8e9f3962-d86c-4e07-b5d3-f31ee5361ca2
animals

# ╔═╡ 430475b3-7db2-4e09-9422-ac1d6ef32ee7
let
	animals
	animals_count[] += 1
end

# ╔═╡ d8613c4f-6936-4fb2-9b9f-acf34377091f
@bind funcs MultiCheckBox([sin, cos, tan]; default=[sin, tan])

# ╔═╡ 92fe974e-cc6c-4387-8c6e-85813b222f25
funcs

# ╔═╡ 1dd8a041-60c0-4b26-8207-3dc9cca0d3eb
[f(0.5) for f in funcs]

# ╔═╡ a8a7e90d-8bbf-4ab6-90a8-24a10885fb0a
@bind animals2 MultiCheckBox(["\"🐰\\\"", "🐱" , "🐵", "🐘", "🦝", "🐿" , "🐝",  "🐪"]; orientation=:column, default=["🐿" , "🐝"])

# ╔═╡ 6123cf6d-fc29-4a8f-a5a1-4366cc6457b6
animals2

# ╔═╡ 60183ad1-4919-4402-83fb-d53b86dda0a6
MultiCheckBox(["🐰 &&\\a \$\$", "🐱" , "🐵", "🐘", "🦝", "🐿" , "🐝",  "🐪"])

# ╔═╡ ad6dcdec-2fc9-45d2-8828-62ac857b4afa
@bind snacks MultiCheckBox(
	["🐱" => "🐝", "🐵" => "🦝", "🐱" => "🐿"];
	default=["🐱", "🐱"]
)

# ╔═╡ abe4c3e0-6e1e-4e26-a4fa-bd60f31c1a4c
snacks

# ╔═╡ Cell order:
# ╟─a8c1e0d2-3604-4e1d-a87c-c8f5b86b79ed
# ╠═8bfaf4c8-557d-433e-a228-aac493746efc
# ╠═8e9f3962-d86c-4e07-b5d3-f31ee5361ca2
# ╠═1c914889-9c98-4aa5-b24f-3027b74feb4a
# ╠═430475b3-7db2-4e09-9422-ac1d6ef32ee7
# ╠═d8613c4f-6936-4fb2-9b9f-acf34377091f
# ╠═92fe974e-cc6c-4387-8c6e-85813b222f25
# ╠═1dd8a041-60c0-4b26-8207-3dc9cca0d3eb
# ╠═a8a7e90d-8bbf-4ab6-90a8-24a10885fb0a
# ╠═6123cf6d-fc29-4a8f-a5a1-4366cc6457b6
# ╠═60183ad1-4919-4402-83fb-d53b86dda0a6
# ╠═abe4c3e0-6e1e-4e26-a4fa-bd60f31c1a4c
# ╠═ad6dcdec-2fc9-45d2-8828-62ac857b4afa
# ╟─c8350f43-0d30-45d0-873b-ff56c5801ac1
# ╟─79c4dc76-efa2-4b7f-99ea-243e3a17f81c
# ╟─499ca710-1a50-4aa1-87d8-d213416e8e30
# ╠═631c14bf-e2d3-4a24-8ddc-095a3dab80ef
# ╠═b65c67ec-b79f-4f0e-85e6-78ff22b279d4
# ╠═430e2c1a-832f-11eb-024a-13e3989fd7c2
# ╠═c38de38d-e900-4309-a9f6-1392af2f245b


text:
using Test
using StaticArrays
using ClimateMachine.VariableTemplates

@testset "VariableTemplates - GPU" begin
    include("test_complex_models_gpu.jl")
end

max_stars_repo_name: mwarusz/CLIMA

max_issues_repo_name: mwarusz/CLIMA

max_forks_repo_name: mwarusz/CLIMA


text:
module DocumenterGitHubActions

# Write your package code here.

end


text:
export
    convert_to_array,
    anynonzero,
    isrealof,
    filter_complex,
    filter_negative,
    ispositive,
    isnegative,
    isnonnegative,
    isnonpositive,
    getcoeff,
    getexponent,
    filterterms,
    findallcoeffs,
    findfirstnonzero,
    findpivotsof,
    nonzeroslicesof,
    dropzeroslices,
    integermultiple

export convert_to_array,
    anynonzero,
    isrealof,
    isnegative,
    filter_complex,
    filter_negative,
    findfirstnonzero,
    findpivotsof,
    nonzeroslicesof,
    dropzeroslices,
    integermultiple

###############################################################################
#                             Types compatibilites                            #
###############################################################################

## Add float arithmetic compatibility with rational type
Base.Float64(x::Nemo.fmpq) = Base.Float64(Rational(x))
Base.Float32(x::Nemo.fmpq) = Base.Float32(Rational(x))
Base.Float16(x::Nemo.fmpq) = Base.Float16(Rational(x))
Base.BigFloat(x::Nemo.fmpq) = Base.BigFloat(Rational(x))

## In case we want to push! a fmpq coeff to a float array
Base.convert(::Type{T}, x::Nemo.fmpq) where {T <: AbstractFloat} = T(Rational(x))
Base.convert(::Type{T}, x::Nemo.fmpz) where {T <: AbstractFloat} = T(Int(x))

## From float to fmpq
(a::Nemo.FlintRationalField)(x::T) where {T <: AbstractFloat} = Nemo.fmpq(rationalize(x))
Nemo.fmpq(x::T) where {T <: AbstractFloat} = Nemo.fmpq(rationalize(x))

## From Nemo integers to Juila rationals
Rational{T}(x::Nemo.fmpz) where {T <: Integer} = Rational{T}(T(x))

## From Polymake matrices to Juila arrays
function convert_to_array(M::Polymake.Matrix, ::T=1) where {T}
    mat = Polymake.@convert_to Matrix{Integer} M
    return T.(transpose(Array(mat)))
end

## From Nemo matrices to Juila arrays
function convert_to_array(M::S, ::T=1) where {S <: MatElem,T}
    return T.(Array(M))
end

###############################################################################
#                         Miscelaneos small functions                         #
###############################################################################

# Nemo.isnonzero(x) = (x) != (zero(x))
# use !iszero instead (it is performance equivalent and more clean)
anynonzero(A) = any((!iszero).(A))

isrealof(rtol) = z -> abs(imag(z)) < rtol #

filter_complex(A, rtol) = real.(filter(isrealof(rtol), A))
filter_negative(A) = filter(isnegative, A)
# filter_positiveorthant(A, rtol) = filternegative(filtercomplex(A, rtol))

ispositive(x) = x>zero(x)
isnegative(x) = x<zero(x)

isnonnegative(x) = !(isnegative(x))
isnonpositive(x) = !(ispositive(x))

function getcoeff(term::Tuple{T,S}) where {T,S}
    return term[1]
end

function getexponent(term::Tuple{T,S}) where {T,S}
    return term[2]
end

function filterterms(predicate, p::MPolyElem)
    return Iterators.filter(predicate, dissect(p))
end

function findallcoeffs(predicate, p::MPolyElem)
    return findall(predicate, collect(coeffs(p)))
end

###############################################################################
#                             Array manipulations                             #
###############################################################################

function findfirstnonzero(A)
    return findfirst((!iszero), A)
end

function findpivotsof(W)
    # return findfirstnonzero.(eachrow(W)) ## it works just for matrices, code below is more generic.
    # return mapslices(findfirstnonzero, W, dims=collect(2:ndims(W)))
    return n -> findfirstnonzero.(eachslice(W, dims=n))
end

function nonzeroslicesof(A)
    # return n->vec(mapslices(anynonzero, M, dims=deleteat!(collect(1:ndims(M)),n)))
    return n -> anynonzero.(eachslice(A, dims=n))
end

function dropzeroslices(A)
    indeces = nonzeroslicesof(A).(1:ndims(A))
    return A[indeces...]
end

"""

    integermultiple(A)

Given an array of rational numbers 'A', returns a multiple 'λA' with integer entries, where 'λ' is the minimal integer with this property.
"""
function integermultiple(A)
    return Int.(abs(lcm(denominator.(A))) .* A)
end


text:
using Proj4
using FactCheck

facts("Proj4") do

println("""
C library version: $(Proj4.version)  [\"$(Proj4._get_release())\"]
geodesic support: $(Proj4.has_geodesic_support)
""")

# Some very basic sanity checking
wgs84 = Projection("+proj=longlat +datum=WGS84 +no_defs")
utm56 = Projection("+proj=utm +zone=56 +south +datum=WGS84 +units=m +no_defs")
nad83 = Projection(Proj4.epsg[4269])

# Reference data computed using GeographicLib's GeoConvert tool
@fact transform(wgs84, utm56, [150.0, -27.0, 0.0]) -->
    roughly([202273.912995055, 7010024.033113679, 0], 1e-6)

dup_wgs84 = Projection("+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs")
@fact transform(wgs84, dup_wgs84, [150 -27 0]) -->
    roughly([150 -27 0], 1e-6)
@fact transform(dup_wgs84, wgs84, [150 -27 0]) -->
    roughly([150 -27 0], 1e-6)

# Taken from https://github.com/proj4js/proj4js/blob/master/test/test.js
sweref99tm = Projection("+proj=utm +zone=33 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs")
rt90 = Projection("+lon_0=15.808277777799999 +lat_0=0.0 +k=1.0 +x_0=1500000.0 +y_0=0.0 +proj=tmerc +ellps=bessel +units=m +towgs84=414.1,41.3,603.1,-0.855,2.141,-7.023,0 +no_defs")
@fact transform(sweref99tm, rt90, [319180 6399862]) -->
    roughly([1271137.927154 6404230.291456], 1e-6)
@fact transform(sweref99tm, rt90, [319480 6397862;
                                   329200 6599800]) -->
    roughly([1271414.272854 6402225.564811;
             1283568.895883 6604160.216834], 1e-6)
@fact transform(sweref99tm, rt90, [319480 6397862;
                                   329200 6599800;
                                   319480 6397862;
                                   329200 6599800]) -->
    roughly([1271414.272854 6402225.564811;
             1283568.895883 6604160.216834;
             1271414.272854 6402225.564811;
             1283568.895883 6604160.216834], 1e-6)
@fact transform(rt90, sweref99tm, [1271414.272854 6402225.564811;
                                   1283568.895883 6604160.216834;
                                   1271414.272854 6402225.564811;
                                   1283568.895883 6604160.216834]) -->
    roughly([319480 6397862;
             329200 6599800;
             319480 6397862;
             329200 6599800], 1e-6)

# Just to verify that we get different Projections from different (projection) definitions
wgs84 = Projection(Proj4.epsg[4326]) # World
svy21 = Projection(Proj4.epsg[3414]) # Singapore
proj = Projection(Proj4.epsg[3906]) # based on the Australian National Spheroid (bessel ellipse)

wgs84_a, wgs84_es = spheroid_params(wgs84)
svy21_a, svy21_es = spheroid_params(svy21)
proj_a, proj_es = spheroid_params(proj)
@fact wgs84_es  --> roughly(svy21_es, 1e-6) # same ellipsoid
@fact svy21_es  --> roughly(0.0066943799901413165, 1e-6)
@fact proj_es   --> roughly(0.006674372231802145, 1e-6)
@fact abs(proj_es - wgs84_es) > 1e-6 --> true
@fact abs(proj_a - wgs84_a) > 500.0 --> true

@fact compare_datums(svy21, wgs84) --> false
@fact compare_datums(svy21, svy21) --> true
@fact is_latlong(wgs84) --> true
@fact is_latlong(utm56) --> false
@fact is_latlong(sweref99tm) --> false
@fact is_latlong(svy21) --> false
@fact is_geocent(Projection(Proj4.epsg[4328])) --> true # WGS 84 (geocentric)
@fact is_geocent(wgs84) --> false # WGS 84 (geocentric)

#=
# TODO(cjf): It's unclear exatly what the below was meant to achieve.
# Commented out for now since the tests are broken.
function test_point(point)
    rev = reverse(vec(point))
    rad = deg2rad(point)
    rev_rad = deg2rad(point)

    @fact lonlat2xy!(rev, wgs84) --> roughly(lonlat2xy(reverse(vec(point)), wgs84), 1e-6) # inplace
    @fact rev --> roughly(rev_rad, 1e-6)
    @fact point --> roughly(rad2deg(rad), 1e-6)
end

for point in ([-73.78, 40.64], # JFK
              [-73.78  40.64],
              [103.99,  1.36], # SIN
              [103.99   1.36])
    test_point(point)
end
=#

p1, p2 = [-73.78, 40.64],[103.99, 30.36] # p1, p2 in degrees
proj1, proj2 = wgs84, proj # chosen to have different ellipses
q1, q2 = lonlat2xy(p1, proj1), lonlat2xy(p2, proj2)

if Proj4.has_geodesic_support
    r2 = transform(proj2, proj1, q2)
    dist_q1r2, azi1, azi2 = geod_inverse(q1, r2, proj1)
    dest, azi = geod_direct(q1, azi1, dist_q1r2, proj1)
    @fact dest --> roughly(r2, 1e-6)
    @fact azi --> roughly(azi2, 1e-6)
    @fact geod_destination(q1, azi1, dist_q1r2, proj1) --> roughly(r2, 1e-6)

    # # It is not necessarily symmetric:
    # # You can go q1(azi1) -[dist1]-> r2,
    # # but not    r2(azi2) -[dist1]-> q1
    # # So the following statements are false in general:
    # dest, azi = geod_direct(r2, azi2, dist1, proj1)
    # @fact dest --> roughly(q1, 1e-6)
    # @fact azi --> roughly(azi1, 1e-6)
    # @fact geod_destination(r2, azi2, dist1, proj1) --> roughly(q1, 1e-6)

    # # To get the reverse azimuth to move from r2 back to q1,
    # # make another call to geod_inverse:
    dist_r2q1, azi1, azi2 = geod_inverse(r2, q1, proj1)
    dest, azi = geod_direct(r2, azi1, dist_r2q1, proj1)
    @fact dest --> roughly(q1, 1e-6)
    @fact azi --> roughly(azi2, 1e-6)
    @fact geod_destination(r2, azi1, dist_r2q1, proj1) --> roughly(q1, 1e-6)

    # The distances from both calls to geod_inverse should be the same, i.e.
    @fact dist_r2q1 --> roughly(dist_q1r2, 1e-6)

    # Doublecheck when we perform the operations in the other projection:
    r1 = transform(proj1, proj2, q1)
    dist_r1q2, azi1, azi2 = geod_inverse(r1, q2, proj2)
    dest, azi = geod_direct(r1, azi1, dist_r1q2, proj2)
    @fact dest --> roughly(q2, 1e-6)
    @fact azi --> roughly(azi2, 1e-6)
    @fact geod_destination(r1, azi1, dist_r1q2, proj2) --> roughly(q2, 1e-6)
    dist_q2r1, azi1, azi2 = geod_inverse(q2, r1, proj1)
    dest, azi = geod_direct(q2, azi1, dist_q2r1, proj1)
    @fact dest --> roughly(r1, 1e-6)
    @fact azi --> roughly(azi2, 1e-6)
    @fact geod_destination(q2, azi1, dist_q2r1, proj1) --> roughly(r1, 1e-6)

    @fact geod_distance(r1, q2, proj2) --> roughly(dist_r1q2, 1e-6)
    @fact geod_distance(q2, r1, proj1) --> roughly(dist_q2r1, 1e-6)
    @fact geod_distance(q1, r2, proj1) --> roughly(dist_q1r2, 1e-6)
    @fact geod_distance(r2, q1, proj1) --> roughly(dist_r2q1, 1e-6)

    # The distances computed in both projections can be significantly different if the ellipsoids are different
    @fact (dist_r1q2 - dist_q1r2) > 5e4 --> true
    @fact (dist_q2r1 - dist_r2q1) > 5e4 --> true

    # So be consistent with your choice of projections when computing distances,
    # and deviate from WGS84 only if you know what you're doing
    @fact geod_distance(r1, q2, proj2) --> roughly(geod_distance(q2, r1, proj2), 1e-6)
    @fact geod_distance(q1, r2, proj1) --> roughly(geod_distance(r2, q1, proj1), 1e-6)
    @fact geod_distance(r1, q2, proj2) --> roughly(geod_distance(q2, r1, proj2), 1e-6)
    @fact geod_distance(q1, r2, proj1) --> roughly(geod_distance(r2, q1, proj1), 1e-6)
end

epsg_error = Int[]
for epsg_code in keys(Proj4.epsg)
    try
        proj = Projection(Proj4.epsg[epsg_code])
        proj_string1 = string(proj)
        proj1 = Projection(proj_string1)
        proj_string2 = string(proj1)
        proj2 = Projection(proj_string2)
        proj_string3 = string(proj2)
        @fact proj_string1 --> proj_string2
        @fact proj_string2 --> proj_string3
        @fact proj_string3 --> proj_string1
    catch
        push!(epsg_error, epsg_code)
    end
end

esri_error = Int[]
for esri_code in keys(Proj4.esri)
    try
        proj_string = Proj4.esri[esri_code]
        proj = Projection(proj_string)
        proj_string1 = string(proj)
        proj1 = Projection(proj_string1)
        proj_string2 = string(proj1)
        proj2 = Projection(proj_string2)
        proj_string3 = string(proj2)
        @fact proj_string1 --> proj_string2
        @fact proj_string2 --> proj_string3
        @fact proj_string3 --> proj_string1
    catch
        push!(esri_error, esri_code)
    end
end

if length(epsg_error) > 0 || length(esri_error) > 0
    errorFraction = (length(epsg_error) + length(esri_error)) /
                    (length(Proj4.epsg) + length(Proj4.esri))
    # Some errors are ok (due to old libproj versions), but a good fraction of
    # the strings should parse - if not something *really* wrong has probably
    # occurred.
    @fact errorFraction --> less_than(0.1)

    println(
    """
    The following projection strings could not be parsed by your version of
    libproj Note that this isn't necessarily a problem, but you won't be able
    to use the projections in question.

    total errors: $(round(100*errorFraction,2))%
    """)
    for epsg_code in sort(epsg_error)
        println("[EPSG:$epsg_code] \"$(Proj4.epsg[epsg_code])\"")
    end
    for esri_code in sort(esri_error)
        println("[ESRI:$esri_code] \"$(Proj4.esri[esri_code])\"")
    end
end

end # facts

FactCheck.exitstatus()


text:
using Snopt
using DiffResults, ForwardDiff

function rosenbrock(x)

    f = (1 - x[1])^2 + 100*(x[2] - x[1]^2)^2

    return f

end

function opt_wrapper(x)
    result = DiffResults.GradientResult(x)
    result = ForwardDiff.gradient!(result, rosenbrock, x)
    f = DiffResults.value(result)
    g = DiffResults.gradient(result)

    fail = false

    c = []
    dcdx = []

    return f, c, g, dcdx, fail
end

x0 = [4.0; 4.0]
lb = [-5.0; -5.0]
ub = [5.0; 5.0]
options = Dict{String, Any}()
options["Derivative option"] = 1
options["Verify level"] = 1
options["Major optimality tolerance"] = 1e-6

xopt, fopt, info = snopt(opt_wrapper, x0, lb, ub, options)


text:
"""
    ProgressiveHedgingSolver(qpsolver::AbstractMathProgSolver; <keyword arguments>)

Return a progressive-hedging algorithm object specified. Supply `qpsolver`, a MathProgBase solver capable of solving quadratic problems.

The following penalty parameter update procedures are available
- `Fixed`:  Fixed penalty (default) ?Fixed for parameter descriptions.
- `adaptive`: Adaptive penalty update ?Adaptive for parameter descriptions.

The following execution policies are available
- `Serial`:  Classical progressive-hedging (default)
- `Synchronous`: Classical progressive-hedging run in parallel ?Synchronous for parameter descriptions.
- `Asynchronous`: Asynchronous progressive-hedging ?Asynchronous for parameter descriptions.

...
# Arguments
- `qpsolver::AbstractMathProgSolver`: MathProgBase solver capable of solving quadratic programs.
- `penalty::AbstractPenalizer = Fixed()`: Specify penalty update procedure (Fixed, Adaptive)
- `execution::AbstractExecuter = Serial`: Specify how algorithm should be executed (Serial, Synchronous, Asynchronous). Distributed variants requires worker cores.
- <keyword arguments>: Algorithm specific parameters, consult individual docstrings (see above list) for list of possible arguments and default values.
...

## Examples

The following solves a stochastic program `sp` created in `StochasticPrograms.jl` using the progressive-hedging algorithm with Ipopt as an `qpsolver`.

```jldoctest
julia> solve(sp,solver=ProgressiveHedgingSolver(IpoptSolver(print_level=0)))
Progressive Hedging Time: 0:00:06 (1315 iterations)
  Objective:  -855.8332803469432
  δ:          9.436947935542464e-7
:Optimal
```
"""
mutable struct ProgressiveHedgingSolver{S <: QPSolver,
                                        E <: Execution,
                                        P <: AbstractPenalizer} <: AbstractStructuredSolver
    qpsolver::S
    execution::E
    penalty::P
    crash::CrashMethod
    parameters::Dict{Symbol,Any}

    function ProgressiveHedgingSolver(qpsolver::QPSolver;
                                      execution::Execution = Serial(),
                                      penalty::AbstractPenalizer = Fixed(),
                                      crash::CrashMethod = Crash.None(), kwargs...)
        S = typeof(qpsolver)
        E = typeof(execution)
        P = typeof(penalty)
        return new{S, E, P}(qpsolver,
                            execution,
                            penalty,
                            crash,
                            Dict{Symbol,Any}(kwargs))
    end
end

function StructuredModel(stochasticprogram::StochasticProgram, solver::ProgressiveHedgingSolver)
    x₀ = solver.crash(stochasticprogram, solver.qpsolver)
    return ProgressiveHedging(stochasticprogram, solver.qpsolver, solver.execution, solver.penalty; solver.parameters...)
end

function add_params!(solver::ProgressiveHedgingSolver; kwargs...)
    push!(solver.parameters, kwargs...)
    for (k,v) in kwargs
        if k ∈ [:qpsolver, :execution, :penalty, :crash]
            setfield!(solver, k, v)
            delete!(solver.parameters, k)
        end
    end
    return nothing
end

function internal_solver(solver::ProgressiveHedgingSolver)
    return get_solver(solver.qpsolver)
end

function optimize_structured!(ph::AbstractProgressiveHedgingSolver)
    return ph()
end

function fill_solution!(stochasticprogram::StochasticProgram, ph::AbstractProgressiveHedgingSolver)
    # First stage
    first_stage = StochasticPrograms.get_stage_one(stochasticprogram)
    nrows, ncols = first_stage_dims(stochasticprogram)
    StochasticPrograms.set_decision!(stochasticprogram, ph.ξ)
    # stochasticprogram.redCosts = try
    #     getreducedcosts(ph.mastersolver.lqmodel)[1:ncols]
    # catch
    #     fill(NaN, ncols)
    # end
    # stochasticprogram.linconstrDuals = try
    #     getconstrduals(ph.mastersolver.lqmodel)[1:nrows]
    # catch
    #     fill(NaN, nrows)
    # end
    # Second stage
    fill_submodels!(ph, scenarioproblems(stochasticprogram))
end

function solverstr(solver::ProgressiveHedgingSolver)
    return "Progressive-hedging solver under $(str(solver.execution)) and $(str(solver.penalty))"
end


text:
module BayesianTools

using Reexport

include("improperdist.jl")
include("Links/Links.jl")
include("ProductDistributions/ProductDistributions.jl")
include("crossmethods.jl")
#export ProductDistribution, link, invlink

end # module

max_stars_repo_name: gragusa/ProductDistributions

max_issues_repo_name: gragusa/ProductDistributions

max_forks_repo_name: gragusa/ProductDistributions


text:
using GBPirate
using Test

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


text:
export NFSPAgent

"""
    NFSPAgent(; rl_agent::Agent, sl_agent::Agent, args...)

Neural Fictitious Self-Play (NFSP) agent implemented in Julia.
See the paper https://arxiv.org/abs/1603.01121 for more details.

# Keyword arguments

- `rl_agent::Agent`, Reinforcement Learning(RL) agent(default `QBasedPolicy` here, use `DQN` for example), which works to search the best response from the self-play process.
- `sl_agent::Agent`, Supervisor Learning(SL) agent(use `BehaviorCloningPolicy` for example), which works to learn the best response from the rl_agent's policy.
- `η`, anticipatory parameter, the probability to use `ϵ-greedy(Q)` policy when training the agent.
- `rng=Random.GLOBAL_RNG`.
- `update_freq::Int`: the frequency of updating the agents' `approximator`.
- `update_step::Int`, count the step.
- `mode::Bool`, used when learning, true as BestResponse(rl_agent's output), false as AveragePolicy(sl_agent's output).
"""
mutable struct NFSPAgent <: AbstractPolicy
    rl_agent::Agent
    sl_agent::Agent
    η
    rng
    update_freq::Int
    update_step::Int
    mode::Bool
end

(π::NFSPAgent)(env::AbstractEnv) = π.sl_agent(env)

RLBase.prob(π::NFSPAgent, env::AbstractEnv, args...) = prob(π.sl_agent.policy, env, args...)

# update env and policy.
function RLBase.update!(π::NFSPAgent, env::AbstractEnv)
    player = current_player(env)
    action = π.mode ? π.rl_agent(env) : π.sl_agent(env)
    π(PRE_ACT_STAGE, env, action)
    env(action)
    π(POST_ACT_STAGE, env, player)
end

(π::NFSPAgent)(stage::PreEpisodeStage, env::AbstractEnv, ::Any) = update!(π.rl_agent.trajectory, π.rl_agent.policy, env, stage)

function (π::NFSPAgent)(stage::PreActStage, env::AbstractEnv, action)
    rl = π.rl_agent
    sl = π.sl_agent

    # update trajectory
    if π.mode
        update!(sl.trajectory, sl.policy, env, stage, action)
        rl(stage, env, action)# also update rl_policy(both learner network and target network).
    else
        update!(rl.trajectory, rl.policy, env, stage, action)
    end

    # update policy
    π.update_step += 1
    if π.update_step % π.update_freq == 0
        if π.mode
            update!(sl.policy, sl.trajectory)
        else
            rl_learn!(rl.policy, rl.trajectory) # only update rl_policy's learner.
            update!(sl.policy, sl.trajectory)
        end
    end
end

function (π::NFSPAgent)(::PostActStage, env::AbstractEnv, player::Any)
    push!(π.rl_agent.trajectory[:reward], reward(env, player))
    push!(π.rl_agent.trajectory[:terminal], is_terminated(env))
end

function (π::NFSPAgent)(::PostEpisodeStage, env::AbstractEnv, player::Any)
    rl = π.rl_agent
    sl = π.sl_agent

    # update trajectory
    # Note that for the `TERMINAL_REWARD` and `SEQUENTIAL` games, some players may not record their real reward and terminated judgment.
    if !rl.trajectory[:terminal][end]
        rl.trajectory[:reward][end] = reward(env, player)
        rl.trajectory[:terminal][end] = is_terminated(env)
    end

    # collect state and dummy action to rl.trajectory
    action = rand(action_space(env, player))
    push!(rl.trajectory[:state], state(env, player))
    push!(rl.trajectory[:action], action)
    if haskey(rl.trajectory, :legal_actions_mask)
        push!(rl.trajectory[:legal_actions_mask], legal_action_space_mask(env, player))
    end

    # update the policy
    π.update_step += 1
    if π.update_step % π.update_freq == 0
        if π.mode
            update!(sl.policy, sl.trajectory)
        else
            rl_learn!(rl.policy, rl.trajectory) # only update rl_policy's learner.
            update!(sl.policy, sl.trajectory)
        end
    end
end

# the supplement function
function rl_learn!(policy::QBasedPolicy, t::AbstractTrajectory)
    # just learn the approximator, not update target_approximator
    learner = policy.learner
    length(t[:terminal]) - learner.sampler.n <= learner.min_replay_history && return

    _, batch = sample(learner.rng, t, learner.sampler)

    if t isa PrioritizedTrajectory
        priorities = update!(learner, batch)
        t[:priority][inds] .= priorities
    else
        update!(learner, batch)
    end
end


text:
module Level
    import ..Log4jl

    bitstype 32 EventLevel
    const levels = Dict{Int32,Symbol}(
        Int32(typemax(Int16)) => :ALL,
        Int32(600) => :TRACE,
        Int32(500) => :DEBUG,
        Int32(400) => :INFO,
        Int32(300) => :WARN,
        Int32(200) => :ERROR,
        Int32(100) => :FATAL,
        Int32(0)   => :OFF
    )
    Base.isless(x::EventLevel, y::EventLevel) = isless(Int32(x), Int32(y))
    Base.convert{T<:Integer}(::Type{T}, x::EventLevel) = convert(T, Intrinsics.box(Int32, x))

    for (v, l) in levels
        eval(:(const $l = Intrinsics.box(EventLevel, $v)))
    end

    function Base.print(io::IO, level::EventLevel)
        vlvl = Int32(level)
        if haskey(levels, vlvl)
            print(io, levels[vlvl])
        end
    end
    Base.show(io::IO,x::EventLevel) = print(io, x, "::EventLevel")

    "Create a custom level and generate convenience functions"
    function add(lvl::Symbol, lvlval::Int32)
        lvlup = symbol(uppercase(string(lvl)))
        lvllw = symbol(lowercase(string(lvl)))
        @assert !isdefined(Level, lvlup) "Level $lvlup is already defined"
        eval(Level, :(const $lvlup = Intrinsics.box(EventLevel, $lvlval)))
        levels[lvlval] = lvlup
        isdefined(Base, lvllw) && eval(Log4jl, :(import Base.$lvllw))
        eval(Log4jl, parse("makemethods(:$lvllw, Level.$lvlup)"))
    end
end
typealias LEVEL Nullable{Level.EventLevel}


text:
module Tst

include("abstract_unittest_functions.jl")

using Test,TimerOutputs,Printf
using HDF5
using JSON
using DataFrames
import Distributions, Random
import Serialization
using MatrixCompletion

# legacy code to be refactored soon!!

# const to = TimerOutput()

#==============================================================================#
#                               TEST OPTIONS                                   #
#==============================================================================#
const TEST_OPTION_PRINT_TIMER    = false
const TEST_OPTION_SMALL_INPUT    = false
const TEST_OPTION_MEDIUM_INPUT   = false
const TEST_OPTION_LARGE_INPUT    = false
const TEST_OPTION_USE_AUTOGRAD   = false
const TEST_OPTION_TEST_R_WRAPPER = false

#==============================================================================#
#                              SUBMODULE FLAGS                                 #
#==============================================================================#
const FLAG_TEST_CONCEPTS                         = false
const FLAG_TEST_SAMPLING                         = false
const FLAG_TEST_MISC                             = false
const FLAG_TEST_RANDOM_OBJECTS                   = false
const FLAG_TEST_DIAGNOSTICS                      = false
const FLAG_TEST_EXPONENTIAL_FAMILY               = false
const FLAG_TEST_INDEXING_TOOLS                   = false
const FLAG_TEST_SPARSE_EIGEN                     = false
const FLAG_TEST_BETTER_MGF                       = false
const FLAG_TEST_ESTIMATOR_MLE                    = false
const FLAG_TEST_ESTIMATOR_MOM                    = false
const FLAG_TEST_MODEL_FITTING                    = false

const FLAG_TEST_LOSS_OPTIMIZER_POISSON           = false
const FLAG_TEST_LOSS_OPTIMIZER_BERNOULLI         = false
const FLAG_TEST_LOSS_OPTIMIZER_GAMMA             = false
const FLAG_TEST_LOSS_OPTIMIZER_GAUSSIAN          = false
const FLAG_TEST_LOSS_OPTIMIZER_NEGATIVE_BINOMIAL = false
const FLAG_TEST_LOSS_OPTIMIZER_MULTINOMIAL       = false
const FLAG_TEST_ALGO_ADMM                        = false
const FLAG_TEST_LIB_MATH                         = false
const FLAG_TEST_PRETTY_PRINTER                   = false
const FLAG_TEST_UTILITY_BATCHUTILS               = false
const FLAG_TEST_SGD_BERNOULLI                    = false
const FLAG_TEST_SGD_GAMMA                        = false

const FLAG_TEST_ALGO_ALM                         = true

# TODO
const FLAG_TEST_ALGO_ADMM_PARALLELL              = false
const FLAG_TEST_ALGO_SVT                         = false
const FLAG_TEST_ALGO_ONEBIT                      = false
const FLAG_TEST_ALGO_OPTSPACE                    = false

#==============================================================================#
#                             SIMULATION FLAGS                                 #
#==============================================================================#

const FLAG_SIMULATION_ADMM_GAMMA                 = false
const FLAG_SIMULATION_ADMM_BERNOULLI             = false
const FLAG_SIMULATION_ADMM_GAUSSIAN              = false
const FLAG_SIMULATION_ADMM_POISSON               = false
const FLAG_SIMULATION_ADMM_GAUSSIAN_BERNOULLI    = false
const FLAG_SIMULATION_ADMM_NEGATIVE_BINOMIAL     = false
const FLAG_SIMULATION_ADMM_MIXED                 = false

#==============================================================================#
#                             VISUALIZATION FLAGS                              #
#==============================================================================#
const FLAG_VISUAL_RANDOM_OBJECTS = false

#==============================================================================#
#                             SIMULATION SCRIPTS                               #
#==============================================================================#
FLAG_SIMULATION_ADMM_GAMMA ?
  include("simulation_runner_gamma.jl")              : @info @sprintf("Skipped: Simulation[vary rank] Gamma")

FLAG_SIMULATION_ADMM_BERNOULLI ?
  include("simulation_runner_bernoulli.jl")          : @info @sprintf("Skipped: Simulation[vary rank] Bernoulli")

FLAG_SIMULATION_ADMM_GAUSSIAN ?
  include("simulation_runner_gaussian.jl")           : @info @sprintf("Skipped: Simulation[vary rank] Gaussian")

FLAG_SIMULATION_ADMM_POISSON ?
  include("simulation_runner_poisson.jl")            : @info @sprintf("Skipped: Simulation[vary rank] Poisson")

FLAG_SIMULATION_ADMM_NEGATIVE_BINOMIAL ?
  include("simulation_runner_negbin.jl")             : @info @sprintf("Skipped: Simulation[vary rank] NegativeBinomial")

FLAG_SIMULATION_ADMM_GAUSSIAN_BERNOULLI ?
  include("simulation_runner_gaussian_bernoulli.jl") : @info @sprintf("Skipped: Simulation Gaussian + Bernoulli")

FLAG_SIMULATION_ADMM_MIXED ?
  include("simulation_runner_mixed.jl")              : @info @sprintf("Skipped: Simulation Mixed")

#==============================================================================#
#                                 TEST SCRIPTS                                 #
#==============================================================================#
FLAG_TEST_MISC ?
  include("test_runner_misc.jl")               : @info @sprintf("Skipped: Miscellaneous Test\n")

FLAG_TEST_RANDOM_OBJECTS ?
  include("test_runner_random_objects.jl")     : @info @sprintf("Skipped: Random Objects Test\n")

FLAG_TEST_SAMPLING ?
  include("test_runner_sampling.jl")           : @info @sprintf("Skipped: Sampling Test\n")

FLAG_TEST_SPARSE_EIGEN ?
  include("test_runner_sparse_eigen.jl")       : @info @sprintf("Skipped: Sparse Eigen Test\n")

FLAG_TEST_INDEXING_TOOLS ?
  include("test_runner_indexing.jl")           : @info @sprintf("Skipped: Indexing Tracker Test\n")

FLAG_TEST_CONCEPTS ?
  include("test_runner_concepts.jl")           : @info @sprintf("Skipped: Concepts Test\n")

FLAG_TEST_DIAGNOSTICS ?
  include("test_runner_diagnostics.jl")        : @info @sprintf("Skipped: Diagnostics Test\n")

FLAG_TEST_EXPONENTIAL_FAMILY ?
  include("test_runner_exponential_family.jl") : @info @sprintf("Skipped: Exponential Family Test\n")

FLAG_TEST_BETTER_MGF ?
  include("test_runner_better_mgf.jl")         : @info @sprintf("Skipped: MGF Test\n")

FLAG_TEST_ESTIMATOR_MLE ?
  include("test_runner_estimator_mle.jl")      : @info @sprintf("Skipped: MLE Test\n")

FLAG_TEST_ESTIMATOR_MOM ?
  include("test_runner_estimator_mom.jl")      : @info @sprintf("Skipped: MOM Test\n")

FLAG_TEST_MODEL_FITTING ?
  include("test_runner_model_fitting.jl")      : @info @sprintf("Skipped: Model Fitting Test\n")

FLAG_TEST_LOSS_OPTIMIZER_POISSON ?
  include("test_runner_poisson_loss.jl")       : @info @sprintf("Skipped: Poisson Loss Test\n")

FLAG_TEST_LOSS_OPTIMIZER_BERNOULLI ?
  include("test_runner_bernoulli_loss.jl")     : @info @sprintf("Skipped: Bernoulli Loss Test\n")

FLAG_TEST_LOSS_OPTIMIZER_GAMMA ?
  include("test_runner_gamma_loss.jl")         : @info @sprintf("Skipped: Gamma Loss Test\n")

FLAG_TEST_LOSS_OPTIMIZER_NEGATIVE_BINOMIAL ?
  include("test_runner_negative_binomial_loss.jl") : @info @sprintf("Skipped: Bernoulli Loss Test\n")

FLAG_TEST_ALGO_ADMM ?
  include("test_runner_admm_small_input.jl")   : @info @sprintf("Skipped: ADMM Small Input Test\n")

FLAG_TEST_LIB_MATH ?
  include("test_runner_lib_math.jl")           : @info @sprintf("Skipped: Math Library Test\n")

FLAG_TEST_PRETTY_PRINTER ?
  include("test_runner_pretty_printer.jl")     : @info @sprintf("Skipped: Pretty Printer\n")

FLAG_TEST_UTILITY_BATCHUTILS ?
  include("test_runner_batch_utils.jl")        : @info @sprintf("Skipped: Utility[Batch Utils]")

FLAG_TEST_SGD_BERNOULLI ?
  include("test_runner_sgd_bernoulli.jl")      : @info @sprintf("Skipped: SGD[Bernoulli]")

FLAG_TEST_SGD_GAMMA ?
    include("test_runner_sgd_gamma.jl")        : @info @sprintf("Skipped: SGD[Gamma]")

FLAG_TEST_ALGO_ALM ?
    include("test_runner_algo_alm.jl")         : @info @sprintf("Skipped: ALGO[ALM]")
#==============================================================================#
#                                VISUAL SCRIPTS                                #
#==============================================================================#

FLAG_VISUAL_RANDOM_OBJECTS ?
  include("visual_random_objects.jl")        : nothing

if TEST_OPTION_PRINT_TIMER
  println()
  println(to)
end

end


text:
@nospecialize
export TagfulPattern, And, Or,
       Literal, Wildcard,
       Deconstrucution, Guard, Effect,
       untagless, TagfulPattern,
       PatternInfo

abstract type TagfulPattern end

struct PatternInfo
    pattern ::TagfulPattern
    typetag :: TypeObject
end

struct And <: TagfulPattern
    ps :: Vector{PatternInfo}
end

struct Or <: TagfulPattern
    ps :: Vector{PatternInfo}
end

struct Literal{T} <: TagfulPattern
    val :: T
end

struct Wildcard <: TagfulPattern
end

struct Deconstrucution <: TagfulPattern
    comp :: PComp
    extract :: Function
    params :: Vector{PatternInfo}
end

struct Guard <: TagfulPattern
    predicate :: Any
end

struct Effect <: TagfulPattern
    perform :: Any
end

@specialize
function _uncurry_call_argtail(f)
    function (_, args...)
        f(args...)
    end
end
@nospecialize

function untagless(points_of_view::Dict{Function, Int})
    myviewpoint = points_of_view[untagless]
    typetag_viewpoint::Int = points_of_view[tag_extract]
    mk_info(all_info)::PatternInfo = PatternInfo(
        all_info[myviewpoint], all_info[typetag_viewpoint]
    )
    ! = mk_info
    function decons(comp::PComp, extract::Function, ps)
        Deconstrucution(comp, extract, PatternInfo[!p for p in ps])
    end
    and(ps::Vector{Vector{Any}}) = And(PatternInfo[!e for e in ps])
    or(ps::Vector{Vector{Any}}) = Or(PatternInfo[!e for e in ps])

    (
        and = and,
        or = or,
        literal = Literal,
        wildcard = Wildcard(),
        decons = decons,
        guard = Guard,
        effect = Effect
    )
end
@specialize


text:
# https://julianlsolvers.github.io/Optim.jl/stable/#examples/generated/ipnewton_basics/

using Optim

fun(x) = (1.0 - x[1])^2 + 100.0 * (x[2] - x[1]^2)^2

function con2_c!(c, x)
    c[1] = x[1]^2 + x[2]^2     ## First constraint
    c[2] = x[2]*sin(x[1])-x[1] ## Second constraint
    return c
end

x0 = [0.25, 0.25]
df = TwiceDifferentiable(fun, x0)

lc = [-Inf, 0.0]; uc = [0.5^2, 0.0]
lx = [-0.5, -0.5]; ux = [0.5, 0.5]
dfc = TwiceDifferentiableConstraints(con2_c!, lx, ux, lc, uc)

res = optimize(df, dfc, x0, IPNewton(), Optim.Options(show_trace=true))
display(res)

println(res.minimum)
println(res.minimizer)

max_stars_repo_name: lanl-ansi/rosetta-opf

max_issues_repo_name: lanl-ansi/rosetta-opf

max_forks_repo_name: lanl-ansi/rosetta-opf


text:
using Glob
using FileIO
using CSV
using DataFrames
using Statistics
using Colors
using Plots
using Plots.PlotMeasures
using LaTeXStrings

#using PackageCompiler
#create_sysimage(:Plots, sysimage_path="sys_plots.so")

pyplot()

include("../config.jl")
include("../utils/utils.jl")

function get_data(datasetname)
    result_path = "../../result/real/$datasetname/"
    horizon_list = glob("$result_path*")

    horizontal_ax = []
    for horizon in horizon_list
        tmp = split(horizon, "/")[end]
        rank_str = split(tmp,"_")[2:end]
        rank = parse.(Int, rank_str)
        #append!(horizontal_ax, rank[1])
        append!(horizontal_ax, prod(rank))
    end

    KLerror = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )
    LSerror = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )
    runtime = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )

    KLerror_std = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )
    LSerror_std = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )
    runtime_std = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )

    for horizon in horizon_list
        for method in methods
            target_path = horizon*"/"*method
            datalists = glob("$target_path*")
            tmp_for_runtime = []
            tmp_for_LSerror = []
            tmp_for_KLerror = []
            for datalist in datalists
                data = load(datalist)
                append!(tmp_for_runtime, data["runtime"])
                append!(tmp_for_LSerror, data["LSerror"])
                append!(tmp_for_KLerror, data["KLerror"])
            end
            append!(LSerror[method], mean( tmp_for_LSerror ))
            append!(KLerror[method], mean( tmp_for_KLerror ))
            append!(runtime[method], mean( tmp_for_runtime ))
            append!(LSerror_std[method], std( tmp_for_LSerror ))
            append!(KLerror_std[method], std( tmp_for_KLerror ))
            append!(runtime_std[method], std( tmp_for_runtime ))
        end
    end
    return horizontal_ax, KLerror, LSerror, runtime, KLerror_std, LSerror_std, runtime_std
end

function plot_error(horizontal_ax, error, error_std, datasetname ;cost="LS", rel=false, xlog=false, std=true)
    ylabel = "$cost error"
    if rel == true
        rel_error = Dict( method for method = zip(methods, [ [] for i=1:length(methods)] ) )
        for method in methods
            rel_error[method] = ( error[method] - error["TLrRR"] ) ./ error[method]
        end
        ylabel = "Relative $cost error"
        error = rel_error
    end

    plt_time = plot()
    if xlog
        plt_time = plot!(plt_time, xaxis=:log)
    end

    title     = "$datasetname"
    xlabel    = "Number of parameters of core tensor"

    #ylim_min = minimum(minimum.(values(error)))
    #ylim_max = maximum(maximum.(values(error)))
    for method in methods
        if std == false
            error_std[method] .= 0.0
        end
        #if method == "TLrRR" && rel == true
        #    continue
        #end
        plot!(plt_time, horizontal_ax[2:end], error[method][2:end],
            title     = title,
            xlabel    = xlabel,
            ylabel    = ylabel,
            grid      = "on",
            legend    = legend,
            legendfont = legendfont,
            size       = img_size,

            #ylim = (ylim_min, ylim_max),

            label       = label_dict[method],
            linestyle   = linetype_dict[method],
            markershapes= markershapes_dict[method],
            linewidth   = linewidth,

            yguidefont  = fnt3,
            xguidefont  = fnt3,
            xtickfont   = fnt3,
            ytickfont   = fnt3,
            #yerror = error_std[method],
            xticks = ([100,1000,10000], [100, 1000, 10000]),
            #yticks = ([25000,35000,45000, 55000], [25000, 35000, 45000, 55000]),

            markersize  = markersize,
            markercolor = :transparent,
            markerstrokewidth = markerstrokewidth,
            markerstrokecolor = colors_dict[method],
            linecolor = colors_dict[method],
            )
    end
    save_path = "../../png/real/$cost"*"error_$datasetname.pdf"
    savefig(plt_time, save_path)
    println("$save_path has been saved")
end

function plot_time(horizontal_ax, runtime, runtime_std, datasetname; xlog, std=true)
    title     = "$datasetname"
    xlabel    = "Number of parameters of core tensor"

    plt_time = plot()
    if xlog
        plt_time = plot!(plt_time, xaxis=:log)
    end

    ylim_min = minimum(minimum.(values(runtime)))
    ylim_max = maximum(maximum.(values(runtime)))
    for method in methods
        if std == false
            runtime_std[method] .= 0.0
        end
        plot!(plt_time, horizontal_ax[2:end], runtime[method][2:end],
            yaxis     = :log,
            title     = title,
            xlabel    = xlabel,
            ylabel    = "Running time(sec.)",
            grid      = "on",
            legend    = legend,
            legendfont = legendfont,
            size       = img_size,
            linewidth  = linewidth,

            ylim = (ylim_min, ylim_max),
            label       = label_dict[method],
            linestyle   = linetype_dict[method],
            markershapes= markershapes_dict[method],

            yguidefont  = fnt3,
            xguidefont  = fnt3,
            xtickfont   = fnt1,
            ytickfont   = fnt1,
            yerror = runtime_std[method],
            xticks = ([100,1000,10000], [100, 1000, 10000]),

            markersize  = markersize,
            markercolor = :transparent,
            markerstrokewidth = markerstrokewidth,
            markerstrokecolor = colors_dict[method],
            linecolor = colors_dict[method],
             )
    end
    save_path = "../../png/real/time_$datasetname.pdf"
    savefig(plt_time, save_path)
    println("$save_path has been saved")
end

function plots_fix_size(datasetname)
    horizontal_ax, KLerror, LSerror, runtime, KLerror_std, LSerror_std, runtime_std = get_data(datasetname)
    println("Data has been collected")
    plot_time(horizontal_ax, runtime, runtime_std, datasetname, xlog=true)
    plot_error(horizontal_ax, LSerror, LSerror_std, datasetname, cost="LS", rel=false, xlog=true, std=false)
    plot_error(horizontal_ax, KLerror, KLerror_std, datasetname, cost="KL", rel=false, xlog=true, std=false)
end

#datasetname = "AttFace"
#plots_fix_size(datasetname)
datasetname = "4DLFD"
plots_fix_size(datasetname)
#datasetname = "TTSB"
#plots_fix_size(datasetname)
#datasetname = "RGBD"
#plots_fix_size(datasetname)
#datasetname = "GTEA"
#plots_fix_size(datasetname)
#datasetname = "fMRI"
#plots_fix_size(datasetname)

max_stars_repo_name: gkazunii/Legendre-tucker-rank-reduction

max_issues_repo_name: gkazunii/Legendre-tucker-rank-reduction

max_forks_repo_name: gkazunii/Legendre-tucker-rank-reduction


text:
using BenchmarkTools
using Dates

"""
    @btimed expression [other parameters...]

An amended version of  BenchmarkTools.@btime. Identical except the return is a tuple of
the result of the `expression` evaluation, the trialmin (of type BenchmarkTools.TrialEstimate)
and the memory allocated (a number of bytes).
"""
macro btimed(args...)
    _, params = BenchmarkTools.prunekwargs(args...)
    bench, trial, result = gensym(), gensym(), gensym()
    trialmin, trialallocs = gensym(), gensym()
    tune_phase = BenchmarkTools.hasevals(params) ? :() : :($BenchmarkTools.tune!($bench))
    return esc(quote
        local $bench = $BenchmarkTools.@benchmarkable $(args...)
        $BenchmarkTools.warmup($bench)
        $tune_phase
        local $trial, $result = $BenchmarkTools.run_result($bench)
        local $trialmin = $BenchmarkTools.minimum($trial)
        local $trialallocs = $BenchmarkTools.allocs($trialmin)
        println("  ",
                $BenchmarkTools.prettytime($BenchmarkTools.time($trialmin)),
                " (", $trialallocs , " allocation",
                $trialallocs == 1 ? "" : "s", ": ",
                $BenchmarkTools.prettymemory($BenchmarkTools.memory($trialmin)), ")")
        $result, $trialmin, $BenchmarkTools.memory($trialmin)
    end)
end

"""
    speedtest(functions, nr::Int, nc::Int)

Prints comparisons of execution speed.
# Arguments
- `functions`:  an array of functions, each an implementation of KendallTau.
- `nr`: number of rows in test matrices.
- `nc`: number of columns in test matrices.

# Example
```
julia>using StatsBase;KendallTau.speedtest([StatsBase.corkendall,KendallTau.corkendall],2000,10)
###################################################################
Executing speedtest 2021-01-19T16:26:29.883
size(matrix1) = (2000, 10)
StatsBase.corkendall(matrix1)
  33.503 ms (451 allocations: 5.54 MiB)
KendallTau.corkendall(matrix1)
  6.172 ms (1918 allocations: 7.82 MiB)
Speed ratio KendallTau.corkendall vs StatsBase.corkendall: 5.428169574078446
Ratio of memory allocated KendallTau.corkendall vs StatsBase.corkendall: 1.4125820189187552
KendallTau.corkendallthreads_v2(matrix1)
  1.878 ms (2234 allocations: 7.86 MiB)
Speed ratio KendallTau.corkendallthreads_v2 vs StatsBase.corkendall: 17.83603066439523
Ratio of memory allocated KendallTau.corkendallthreads_v2 vs StatsBase.corkendall: 1.4198018874681153
Results from all 3 functions identical? true
--------------------------------------------------
```
"""
function speedtest(functions, nr::Int, nc::Int)

    rng = MersenneTwister(1)# make the contents of matrix1 etc. deterministic so successive calls with fixed nc & nr are fully comparable
    results = Array{Any}(undef, length(functions))
    times = Array{Float64}(undef, length(functions))
    allocations = Array{Float64}(undef, length(functions))
    matrix1 = randn(rng, Float64, nr, nc)
    matrix2 = randn(rng, Float64, nr, nc)
    vector1 = randn(rng, nr)
    vector2 = randn(rng, nr)
    manyrepeats1 = rand(rng, 1:2, nr)
    manyrepeats2 = rand(rng, 1:2, nr)

    fname = f -> string(Base.parentmodule(f)) * "." * string(f)

    println("#"^67)
    println("Executing speedtest $(now())")

    for k = 1:5
        println("-"^50)
        if k == 1
            @show(size(matrix1))
        elseif k == 2
            @show(size(matrix1))
            @show(size(matrix2))
        elseif k == 3
            @show(size(vector1))
            @show(size(matrix1))
        elseif k == 4
            @show(size(matrix1))
            @show(size(vector1))
        elseif k == 5
            @show(size(vector1))
            @show(size(vector2))
        end
        i = 0
        for fn in functions
            i += 1
            if k == 1
                println("$(fname(fn))(matrix1)")
                tmp =  @btimed $fn($matrix1)
            elseif k == 2
                println("$(fname(fn))(matrix1,matrix2)")
                tmp =  @btimed $fn($matrix1, $matrix2)
            elseif k == 3
                println("$(fname(fn))(vector1,matrix1)")
                tmp =  @btimed $fn($vector1, $matrix1)
            elseif k == 4
                println("$(fname(fn))(matrix1,vector1)")
                tmp =  @btimed $fn($matrix1, $vector1)
            elseif k == 5
                println("$(fname(fn))(vector1,vector2)")
                tmp =  @btimed $fn($vector1, $vector2)
            end
            results[i], times[i], allocations[i] = tmp[1], tmp[2].time, tmp[3]
            if i > 1
                println("Speed ratio $(fname(functions[i])) vs $(fname(functions[1])): $(times[1] / times[i])")
                println("Ratio of memory allocated $(fname(functions[i])) vs $(fname(functions[1])): $(allocations[i] / allocations[1])")
            end
        end
        if length(functions) > 1
          nfns = length(functions)
            resultsidentical = all(myapprox.(results[2:end], results[1:(end - 1)], 1e-14))
            if !resultsidentical
                @warn "Results from $(nfns ==2 ? "both" : "all $nfns") functions identical? $resultsidentical"
            else
                println("Results from $(nfns ==2 ? "both" : "all $nfns") functions identical? $resultsidentical")
            end
        end
    end
    println("#"^67)
end

# Test for equality with absolute tolerance of `abstol` and NaN being equal to NaN (different to Julia's isapprox)
function myapprox(x::AbstractArray, y::AbstractArray, abstol::Float64)
    if size(x) ≠ size(y)
        return(false)
    else
        return(all(myapprox.(x, y, abstol)))
    end
end

function myapprox(x::Float64, y::Float64, abstol::Float64)
    if isnan(x) && isnan(y)
        return(true)
    elseif isnan(x)
        return(false)
    elseif isnan(y)
        return(false)
    else
        return(abs(x - y) <= abstol)
    end
end

# Code to investigate performance impact of the presence of missings in the arguments passed to corkendall
function sprinklemissings(x, proportionmissing, rng=MersenneTwister())
    randoms = rand(rng, size(x)...)
    ifelse.(randoms .< proportionmissing, missing, x)
end

function impactofmissings(nr::Int, nc::Int, proportionmissing::Float64=0.1)

    fn1 = KendallTau.corkendall
    fn2 = KendallTau.corkendall# BaseStats.corkendall

    rng = MersenneTwister(1)# make the contents of matrix1 etc. deterministic so successive calls with fixed nc & nr are fully comparable
    results = Array{Any}(undef, 2)
    times = Array{Float64}(undef, 2)
    allocations = Array{Float64}(undef, 2)
    matrix1 = randn(rng, Float64, nr, nc)
    matrix2 = randn(rng, Float64, nr, nc)
    vector1 = randn(rng, nr)
    vector2 = randn(rng, nr)

    fname = f -> string(Base.parentmodule(f)) * "." * string(f)

    println("#"^67)
    println("Executing impactofmissings $(now())")

    for k = 1:5
        println("-"^50)
        if k == 1
            @show(size(matrix1))
        elseif k == 2
            @show(size(matrix1))
            @show(size(matrix2))
        elseif k == 3
            @show(size(vector1))
            @show(size(matrix1))
        elseif k == 4
            @show(size(matrix1))
            @show(size(vector1))
        elseif k == 5
            @show(size(vector1))
            @show(size(vector2))
        end
        i = 0
        functions = [fn1,fn2]
        for fn in functions
            i += 1
            message = " no missings in argument(s)"
            if i == 2
                message = " argument(s) amended to contain $(100proportionmissing)% missings"
                matrix1 = sprinklemissings(matrix1, proportionmissing, rng)
                matrix2 = sprinklemissings(matrix2, proportionmissing, rng)
                vector1 = sprinklemissings(vector1, proportionmissing, rng)
                vector2 = sprinklemissings(vector2, proportionmissing, rng)
            end

            if k == 1
                println("$(fname(fn))(matrix1)$message")
                tmp =  @btimed $fn($matrix1)
            elseif k == 2
                println("$(fname(fn))(matrix1,matrix2)$message")
                tmp =  @btimed $fn($matrix1, $matrix2)
            elseif k == 3
                println("$(fname(fn))(vector1,matrix1)$message")
                tmp =  @btimed $fn($vector1, $matrix1)
            elseif k == 4
                println("$(fname(fn))(matrix1,vector1)$message")
                tmp =  @btimed $fn($matrix1, $vector1)
            elseif k == 5
                println("$(fname(fn))(vector1,vector2)$message")
                tmp =  @btimed $fn($vector1, $vector2)
            end
            results[i], times[i], allocations[i] = tmp[1], tmp[2].time, tmp[3]
            if i > 1
                println("Speed ratio $(fname(functions[i])) ($(100proportionmissing)% missings) vs $(fname(functions[1])) (no missings): $(times[1] / times[i])")
                println("Ratio of memory allocated $(fname(functions[i])) ($(100proportionmissing)% missings) vs $(fname(functions[1])) (no missings): $(allocations[i] / allocations[1])")
            end
        end
    end
    println("#"^67)
end

using Missings

function sm1(x, y)
    mx, my = Missings.skipmissings(x, y)
    collect(mx), collect(my)
end

function testmissings()
    x = [missing;1:1000]
    y = [1:1000;missing]
    @benchmark skipmissingpairs($x, $y)
end

#= test different ways of "skipping missing pairs".
julia> KendallTau.test_skipmissings(10000)
  46.700 μs (7 allocations: 181.58 KiB)
  151.800 μs (46 allocations: 514.88 KiB)
  25.400 μs (4 allocations: 156.41 KiB) =#
function test_skipmissings(n=10000)

    x = [missing;1:n]
    y = [1:n;missing]

    # simplest approach I could think of
    @btime begin
        keep = .!(ismissing.($x) .| ismissing.($y))
        x2 = $x[keep]
        y2 = $y[keep]
    end

    # using Missings.skipmissings
    @btime begin
        itrx, itry = Missings.skipmissings($x, $y)
        # I think I may be misusing Missings.skipmissings by calling collect here
        x3 = collect(itrx)
        y3 = collect(itry)
    end

    # use KendallTau.skipmissingpairs
    @btime x4, y4 = KendallTau.skipmissingpairs($x, $y)

    nothing
end


text:
incrmark(x::Float64)   = setmark(x, Marker1)
incrmark(x::MarkableF640) = setmark(x, Marker1)
incrmark(x::MarkableF641) = setmark(x, Marker2)
incrmark(x::MarkableF642) = setmark(x, Marker3)
incrmark(x::MarkableF643) = x

decrmark(x::Float64)  = x
decrmark(x::MarkableF640) = x
decrmark(x::MarkableF641) = setmark(x, Marker0)
decrmark(x::MarkableF642) = setmark(x, Marker1)
decrmark(x::MarkableF643) = setmark(x, Marker2)

incrmark(x::Float32)   = setmark(x, Marker1)
incrmark(x::MarkableF320) = setmark(x, Marker1)
incrmark(x::MarkableF321) = setmark(x, Marker2)
incrmark(x::MarkableF322) = setmark(x, Marker3)
incrmark(x::MarkableF323) = x

decrmark(x::Float32)  = x
decrmark(x::MarkableF320) = x
decrmark(x::MarkableF321) = setmark(x, Marker0)
decrmark(x::MarkableF322) = setmark(x, Marker1)
decrmark(x::MarkableF323) = setmark(x, Marker2)

incrmark(x::Float16)   = setmark(x, Marker1)
incrmark(x::MarkableF160) = setmark(x, Marker1)
incrmark(x::MarkableF161) = setmark(x, Marker2)
incrmark(x::MarkableF162) = setmark(x, Marker3)
incrmark(x::MarkableF163) = x

decrmark(x::Float16)  = x
decrmark(x::MarkableF160) = x
decrmark(x::MarkableF161) = setmark(x, Marker0)
decrmark(x::MarkableF162) = setmark(x, Marker1)
decrmark(x::MarkableF163) = setmark(x, Marker2)


text:
"""
    Direct

This module exports all the direct interface functions. This is a low-level interface to the
Dobot that directly exposes all Dobot functionality.
"""
module Direct

using ..DobotMagician: DobotMagician

# Automatically export all Commands
for name in names(DobotMagician; all=true)
    if getproperty(DobotMagician, name) isa DobotMagician.Command
        eval(:(const $name = DobotMagician.$name))
        eval(:(export $name))
    end
end

end


text:
using Printf
const N2048 = 4
const NS2048 = N2048^2
const status = zeros(Int, 4, 4)

const greet_string = """
This is a 2048 game;
use 'w'(up), 's'(down), 'a'(left), 'd'(right) to play.
Ready? [Y/N]"""

include("console.jl")

function start_2048()
    print(greet_string)
    choose = getch()
    if choose != 'Y' % Cint && choose != 'y' % Cint
        return
    end
    status[mod1(rand(Int), NS2048)] = 2
    clean_console()
    print_frame()
    while step2048() end
end

function print_frame()
    set_cursor_position(0, 0)
    line = repeat('-', 30)
    println(line)
    for i = 1:N2048
        for j = 1:N2048
            print("| ")
            if status[i, j] != 0
                @printf("%4d", status[i,j])
            else
                print("    ")
            end
            print(" ")
        end
        println(" |")
        println(line)
    end
end

function shift2048(vec::Vector{Int})::Vector{Int}
    N = size(vec, 1)
    temp = zeros(Int, N)
    shift = zeros(Int, N)
    mask = vec .!= 0
    read_idx = count(mask)
    temp[1:read_idx] .= vec[mask]
    write_idx = 1
    i = 1
    while i <= read_idx
        if i < N && temp[i] == temp[i + 1]
            shift[write_idx] = temp[i] * 2
            write_idx += 1
            i += 2
        else
            shift[write_idx] = temp[i]
            write_idx += 1
            i += 1
        end
    end
    return shift
end

function change_status(direction::Symbol)
    if direction == :up
        for i = 1:N2048
            status[:, i] .= shift2048(status[:, i])
        end
    elseif direction == :down
        for i = 1:N2048
            status[end:-1:1, i] .= shift2048(status[end:-1:1, i])
        end
    elseif direction == :left
        for i = 1:N2048
            status[i, :] .= shift2048(status[i, :])
        end
    elseif direction == :right
        for i = 1:N2048
            status[i, end:-1:1] .= shift2048(status[i, end:-1:1])
        end
    end
    # println("after shift")
    # print_frame()
end

function add_new()
    zero_counts = count(status .!= 0)
    zero_counts == 0 && return
    rand_num = mod(rand(Int), zero_counts)
    for i = 1:NS2048
        zero_counts -= status[i] == 0
        if zero_counts == rand_num
            r = rand(Float64)
            max_s = maximum(status)
            tab = Dict(
                    64 => floor(Int, r * 8.),
                    128 => floor(Int, r * 5.),
                    256 => floor(Int, r * 3.),
                    512 => floor(Int, r * 2.),
                    1024 => floor(Int, r * 1.5)
                )
            p = one(Int)
            if max_s in keys(tab)
                p = tab[max_s]
            end
            if p == 0
                status[i] = 4
            else
                status[i] = 2
            end
            return
        end
    end
end

function step2048()
    ch = getch()
    if ch in ['w' % Int, 72, 65]
        change_status(:up)
    elseif ch in ['s' % Int, 80, 66]
        change_status(:down)
    elseif ch in ['a' % Int, 75, 68]
        change_status(:left)
    elseif ch in ['d' % Int, 77, 67]
        change_status(:right)
    elseif ch in ['q' % Int, 'Q' % Int]
        print("You quit!")
        return false
    else
        return true
    end
    check = check2048()
    if check == 1
        print_frame()
        print("You win!!!")
        return false
    elseif check == -1
        print_frame()
        print("You fail!")
        return false
    end
    add_new()
    print_frame()
    return true
end

function check2048()::Int
    if any(status .== 2048)
        return 1
    elseif any(status .== 0)
        return 0
    elseif any(status[:, 1:end - 1] .== status[:, 2:end])
        return 0
    elseif any(status[1:end - 1, :] .== status[2:end, :])
        return 0
    else
        return -1
    end
end

start_2048()