text:
module EigenSysTests
using Compat

# test eigenSys
include(joinpath(dirname(@__FILE__), "..", "src", "AndersonMoore.jl"))
using .AndersonMoore

#tweaked= False
# test eigenSys! firmvalue example
function firmvalue()::Bool

neq=4::Int64;nlag=1::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=2

aa=[1.1 -0.4 0.;
0. 0.4 0.;
0. 0. 1.]::Array{Float64,2}

ww=[0.868243 0. 0.;
-0.496139 0. 1.;
0. 1. 0.]::Array{Float64,2}

rts=[1.1+0.0im;
     1.0+0.0im;
     0.4+0.0im]::Array{Complex{Float64},1}
lgroots=1

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rts; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=1e-15::Float64) &&
isapprox(rtsJulia,rts,rtol=0.0::Float64,atol=0.1e-10::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! firmvalue3Leads2Lags example
function firmvalue3Leads2Lags()::Bool

neq=4::Int64;nlag=2::Int64;nlead=3::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=4

aa=[0. 1. 0. 0. 0. 0.;
0. 0. 0. 1. 0. 0.;
0. 0. 0. 0. 1. 0.;
1.331 0. -1. 0. 0. 0.;
0. 0. 0.4 0. 0. 0.;
0. 0. 0. 0. 0. 1.]::Array{Float64,2}

ww=[0.455955 0.556804 0.556804 0. 0.0000000000000000840192 0.0000000000000000132968;
0.414505 -0.691462 0.185277 0. 0.0000000000000000700403 -0.0000000000000000285265;
-0.511734 -0.42174 -0.255383 0. -0.534522 0.534522;
0.376823 0.168433 -0.628602 0. -0.0000000000000000500741 0.000000000000000029589;
-0.465213 0.392762 -0.21595 0. 0.845154 0.845154;
0. 0. 0. 1. 0. 0.]::Array{Float64,2}

rts=[1.1+0.0im;
    -0.55+0.952628im;
    -0.55-0.952628im;
     1.0+0.0im;
    -0.632456+0.0im;
     0.632456+0.0im]::Array{Complex{Float64},1}

lgroots=3

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! example7 example
function example7()::Bool

neq=4::Int64;nlag=1::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=3

aa=[1.30303 -1.0101 1.;
-0.30303 1.0101 0.;
-0.333333 1.11111 0.66]::Array{Float64,2}

ww=[-0.0975567 -0.665933 -0.136197;
0.663306 0.663306 -0.897059;
0.391168 -0.716975 0.420399]::Array{Float64,2}

rts=[1.31855+0.49536im;
1.31855-0.49536im;
0.336028+0.0im]::Array{Complex{Float64},1}

lgroots=2

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! oneEquationNoLead example
function oneEquationNoLead()::Bool

neq=1::Int64;nlag=1::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=1

aa=hcat([-0.666667])::Array{Float64,2}

ww=hcat([1.])::Array{Float64,2}

rts=[-0.666667+0.0im]::Array{Complex{Float64},1}

lgroots=0

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! reliablePaperExmpl example
function reliablePaperExmpl()::Bool

neq=5::Int64;nlag=1::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=3

aa=[0. 0. 0. 0.;
-1.11111 1.11111 -0.444444 0.;
-1.55556 1.55556 2.37778 -2.;
0. 0. 2. -1.]::Array{Float64,2}

ww=[0.548572 0.548572 0.604781 1.;
-0.548572 -0.548572 -0.604781 0.;
-0.475037 0.508387 0.291542 0.;
0.544611 -0.124726 -0.428349 0.]::Array{Float64,2}

rts=[1.06383+1.39432im;
1.06383-1.39432im;
0.361236+0.0im;
0.0+0.0im]::Array{Complex{Float64},1}

lgroots=2

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! athan example
function athan()::Bool

neq=9::Int64;nlag=1::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=2

aa=[0. 1. 0. 0. 0. 0. 0.;
-1. 2. -0.36 0. 0. -2. 0.;
0. 0.293 0.764 -0.293 0. 0. 0.;
-2. 4.293 0.044 -0.293 -0.911 -4. 0.;
0. 0. 0. 0. 0.911 0. 0.;
0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 1.]::Array{Float64,2}

ww=[-0.381771 0. 0. 0.330041 0.330041 0.000000000000000161051 0.;
0.447207 0. 0. -0.454296 0.269733 -0.753484 0.;
-0.250613 0. 0. -0.147243 0.605615 -0.376742 0.;
0.0434036 0. 0. 0.102035 -0.186692 0.376742 0.;
-0.081062 0. 1. -0.015374 -0.444561 0.376742 0.;
-0.763541 0. 0. 0.660083 0.660083 0.0803643 1.;
0. 1. 0. 0. 0. 0. 0.]::Array{Float64,2}

rts=[1.39878+0.0im;
1.0+0.0im;
0.911+0.0im;
0.536108+0.222056im;
0.536108-0.222056im;
0.000000000000000541234+0.0im;
0.0+0.0im]::Array{Complex{Float64},1}

lgroots=1

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

#tweaked= False
# test eigenSys! habitmod example
function habitmod()::Bool

neq=12::Int64;nlag=4::Int64;nlead=1::Int64
qRows=(neq*nlead)::Int64;qCols=(neq*(nlag+nlead))::Int64

uprbnd=1.
rowsLeft=6

aa=[0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0.;
-0.155117 -0.0214875 -0.0117771 0.00589435 0.185738 0.00238782 0.0597264 -0.0268822 0.00412584 0.0149078 -0.0493113 0.0238367 -0.872076 0.88879 0.89334 0.100478 -0.103176 -3.54054 36.5013 1.42383 0.744321 -0.00209401 -0.49621 0.114628;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.88879 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.;
-0.52802 -0.11279 -0.031833 -0.024734 0.71526 -0.042076 0.22829 0.012267 -0.22215 0.14172 -0.26679 0.11959 0.42308 0. 0.61076 0.10946 -0.17172 0. 0. 0. 0. 0. 0. 0.00011277;
-0.46436 -0.16131 -0.02744 -0.17087 1.3425 -0.16794 0.35818 0.09837 -1.06382 0.02255 -0.51594 -0.08949 0.46672 0. 0.11369 1.10657 0.26021 0. 0. 0. 0. 0. 0. 0.00141;
0.23292 -0.08512 0.04213 -0.14149 -0.03313 -0.15223 -0.01572 0.34762 -0.70552 0.24697 -0.21888 0.07805 0.56316 0. 0.08268 0.16054 0.68637 0. 0. 0. 0. 0. 0. 0.00348;
0.005823 -0.002128 0.00105325 -0.00353725 -0.00082825 -0.00380575 -0.000393 0.0086905 -0.017638 0.00617425 -0.005472 0.00195125 0.014079 0. 0.002067 -0.0209865 0.0171593 1.025 0. 0. 0. 0. 0. 0.000087;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.99853 0. 0. 0. 0. 0. 0.00147 0. 0. 0. 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. -681.515 0. 0. 0. 0. 0. 17667.9 680.953 0. 0. 0. 0.;
0.530258 0.113268 0.0319679 0.0248388 -0.718291 0.0422543 -0.229257 -0.012319 0.223091 -0.142321 0.267921 -0.120097 -0.424873 0. 0.39089 -0.109924 0.172448 0. 0. 0. 1.00424 0. 0. -0.000113248;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 684.403 0. 0. 0. 0. 0. -17742.8 -682.835 0. 1.00424 0. 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. -1.00276 0. 0. 0. 0. 0. 1.00276 0. 0. 0. 1.00424 0.;
0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1.]::Array{Float64,2}

ww=[0.00000847733 -0.0489971 0.00329698 0.00364315 -0.00296384 0. 0.00289826 0.108025 0. -0.0438951 0.017043 0.0201472 -0.0791448 -0.0676478 0.204738 0.00741613 0.0457941 -0.0564212 -0.00691027 -0.00120837 -0.000000000104862 0.000000000000784312 -0.000000000000767084 -0.000000000106419;
0.00000117376 -0.0112006 0.000706552 0.000780736 0.0000660393 0. -0.0176289 -0.0074902 0. -0.0033955 0.00895106 0.00876825 -0.0129104 -0.0103006 0.028083 0.0304513 0.00587952 -0.0187752 -0.000579761 0.00136933 0.0000000000000976943 0.0000000000000066158 -0.00000000000000641318 0.0000000000000846545;
0.00000064385 -0.00419088 0.000225281 0.000248934 -0.00028174 0. 0.000265449 0.0158931 0. -0.00228812 0.000766883 0.00000997697 -0.00614689 -0.00623485 0.0148889 -0.00363062 0.00467283 0.00161117 -0.00049169 -0.000532766 -0.0000000000000252739 0.00000000000000444306 -0.00000000000000435659 -0.0000000000000340815;
-0.000000322128 -0.0119854 0.00032908 0.000363632 0.000922354 0. -0.0441305 -0.010256 0. 0.0149254 0.0121877 0.00423626 -0.00265616 -0.00671859 -0.0126441 0.0452943 0.00474333 0.0127896 0.0014004 0.00101101 -0.00000000000000494435 -0.00000000000000500096 0.00000000000000488038 0.00000000000000494495;
-0.0000101418 0.0895547 -0.00245965 -0.0027179 -0.0024055 0. 0.200931 -0.107545 0. -0.0900449 -0.0628586 0.0156808 0.254502 0.0432761 -0.397022 0.133855 -0.0655314 0.115807 -0.0154165 -0.0167395 0.000000719464 -0.00000074944 0.000000749437 -0.00000078144;
-0.000000129057 -0.0207267 0.00108315 0.00119687 0.000894541 0. -0.0567803 -0.0272877 0. 0.00609649 0.0186839 0.0242131 0.0165386 -0.0222346 -0.0754314 0.0606591 0.0329998 0.047426 -0.00405635 0.00410942 -0.000000000682305 0.00000000064726 -0.000000000647247 0.000000000612164;
-0.00000326377 0.0311734 -0.00146039 -0.00161372 0.000153904 0. 0.0443083 -0.039071 0. -0.00331164 -0.0177129 -0.00839306 0.0485082 0.0314834 -0.0792245 -0.0148154 -0.0394007 0.0012426 0.00120493 -0.00215158 0.000000000222065 -0.000000000210693 0.000000000210689 -0.00000000019929;
0.00000147074 -0.0160863 0.000408635 0.000451539 -0.000584756 0. -0.0213281 0.0800915 0. 0.0159632 0.00532987 -0.0000349571 -0.000670496 -0.019063 -0.0240676 -0.121402 -0.0196196 0.0286035 0.000173991 0.00151099 0.0000000000000245739 -0.00000000000000517572 0.00000000000000502549 0.0000000000000347826;
-0.000000239717 -0.000792326 0.0000379058 0.0000418857 0.0022348 0. -0.0287681 -0.110363 0. -0.0442625 -0.0259932 0.00592813 -0.180197 0.323735 0.106488 -0.446833 -0.091663 -0.259238 0.00224181 -0.0837174 0.00034508 0.0102842 0.0110113 -0.000383324;
-0.000000811654 -0.0311743 0.000721688 0.000797462 0.00092268 0. -0.0906195 0.0653455 0. 0.0382729 0.0108365 0.0250779 0.0303865 0.0158253 -0.0580298 -0.193109 -0.139833 0.0894769 -0.0202529 -0.00355066 0.000000123028 -0.00000931441 -0.00000906483 -0.000000126568;
0.00000268973 -0.01287 0.00049266 0.000544387 0.00083439 0. -0.0300244 -0.0284717 0. 0.00447302 0.0108379 -0.00467599 -0.0998469 0.0278434 0.146725 -0.00860273 0.0905782 -0.0943029 0.00879296 -0.00123068 -0.0000000398327 0.00000303173 0.00000295089 0.000000041007;
-0.00000129793 -0.0312023 0.000213646 0.000236078 0.000286668 0. -0.0907005 0.151454 0. 0.0594932 0.0159445 -0.0111903 0.0325489 -0.01021 0.153668 0.0785177 -0.0105616 0.00218828 -0.000218761 0.00114412 -0.000000000222072 0.000000000210686 -0.000000000210682 0.000000000199297;
-0.0370804 0.0246579 -0.00000000000000307301 0.00000000000000369367 0.00016314 0. 0.00390901 0.042236 0. -0.0708795 -0.0430473 -0.17663 -0.121307 -0.402926 0.168196 0.378419 0.222077 0.350449 0.162999 0.244079 -0.252057 -0.252092 0.252023 0.252057;
-0.0000485593 0.0958625 -0.0000000000000235343 0.0000000000000271351 0.0012582 0. -0.369591 0.161292 1. 0.361963 0.329714 0.229698 0.23356 0.209559 -0.209672 -0.194319 -0.231622 -0.219961 -0.26965 -0.250136 0.252057 0.252052 -0.252063 -0.252058;
-0.0000488121 0.127132 0.489846 0.541277 -0.5851 0. -0.434604 -0.00160089 0. 0.346801 0.203295 0.183594 0.193904 0.218294 -0.0525509 -0.0682873 -0.110412 -0.370963 -0.229611 -0.256183 0.252057 0.252066 -0.252049 -0.252057;
-0.00000549337 0.12544 -0.00218918 -0.00241903 -0.00344982 0. 0.253458 -0.250185 0. -0.0925811 -0.0341928 0.0312099 0.05905 -0.0406851 -0.096241 0.0167473 -0.061666 0.0923319 -0.000184476 0.00342656 -0.0000000000000157287 -0.0000000000000159155 0.0000000000000155232 0.0000000000000157294;
0.00000563828 -0.0259403 0.000991241 0.00109532 -0.0042903 0. 0.00882998 0.295823 0. 0.0250455 0.00230957 -0.0163765 -0.0328191 0.0263861 -0.10227 0.0107738 0.0347595 -0.0520895 -0.010955 -0.0036157 0.0000000398328 -0.00000303173 -0.00000295089 -0.000000041007;
0.000193477 -0.944771 0.00000000000192237 -0.00000000000102525 0.0275468 0. 0.764362 0.764362 0. -0.808862 -0.888581 -0.888581 -0.804108 -0.804108 0.733356 0.733356 0.807099 0.807099 0.872336 0.872336 -0.870652 -0.870652 0.870653 0.870653;
0.998571 -0.109495 -0.509098 -0.453563 0.393258 0. 0.132405 0.104521 0. -0.117758 -0.129437 -0.116763 -0.10817 -0.0916235 0.0952083 0.0831025 0.102589 0.0954434 0.117021 0.112195 -0.000346297 -0.0102456 -0.0109752 0.000384617;
0.0384857 -0.000509914 -0.000741474 0.0035434 0.00166552 0. 0.000553558 0.000363793 0. -0.000366298 -0.000468001 -0.000164372 -0.000233277 0.000362695 0.0000913327 -0.000352556 0.0000236021 -0.000248118 0.000192018 0.0000171507 -0.000000000984304 0.0000000749258 0.0000000729283 0.00000000101334;
-0.0000406731 0.162178 0.49167 0.543293 -0.589634 0. -0.198838 -0.160111 0. 0.182241 0.197634 0.189195 0.172734 0.170126 -0.157146 -0.155475 -0.173005 -0.171768 -0.188855 -0.187001 0.18682 0.186819 -0.18682 -0.18682;
0.000000114426 -0.000456259 -0.000738341 0.00352842 0.00165883 0. 0.000559392 0.000450443 0. -0.000512702 -0.000556006 -0.000532264 -0.000485953 -0.000478617 0.0004421 0.000437399 0.000486716 0.000483238 0.000531309 0.000526092 -0.000525582 -0.000525581 0.000525584 0.000525583;
0.0000271152 -0.108118 -0.50902 -0.453954 0.393087 0. 0.132557 0.10674 0. -0.121493 -0.131755 -0.126129 -0.115155 -0.113416 0.104763 0.103649 0.115336 0.114511 0.125903 0.124666 -0.124546 -0.124545 0.124546 0.124546;
-0.00000626373 0.0240777 -0.0000183887 -0.0000203194 0.000519754 1. -0.0390233 -0.026122 0. 0.0288168 0.031094 0.0292324 0.0267982 0.0262171 -0.0240533 -0.0241842 -0.0268358 -0.0264974 -0.0293165 -0.028979 0.0289614 0.0289613 -0.0289615 -0.0289614]::Array{Float64,2}

rts=[679.579+0.0im;
1.11741+0.0im;
1.00424+0.0im;
1.00424+0.0im;
1.00403+0.0im;
1.0+0.0im;
0.918128+0.0887657im;
0.918128-0.0887657im;
0.88879+0.0im;
0.714747+0.0im;
0.587349+0.26613im;
0.587349-0.26613im;
-0.129115+0.619993im;
-0.129115-0.619993im;
-0.294489+0.507462im;
-0.294489-0.507462im;
-0.230764+0.281568im;
-0.230764-0.281568im;
0.198899+0.16454im;
0.198899-0.16454im;
-0.0000704459+0.0000723217im;
-0.0000704459-0.0000723217im;
0.0000704459+0.0000685148im;
0.0000704459-0.0000685148im]::Array{Complex{Float64},1}

lgroots=5

(wwJulia,rtsJulia,lgrootsJulia)=eigenSys!(aa,uprbnd,rowsLeft)

    # rounding to 16 digit places
wwJulia = Compat.round.(wwJulia; digits=16, base=10)
ww = Compat.round.(wwJulia; digits=16, base=10)
rtsJulia = Compat.round.(rtsJulia; digits=16, base=10)
rts = Compat.round.(rtsJulia; digits=16, base=10)

isapprox(wwJulia,ww,rtol=0.::Float64,atol=0.1e-16::Float64)&&
isapprox(rtsJulia,rts,rtol=0.::Float64,atol=0.1e-16::Float64)&&
lgrootsJulia==lgroots
end;

end


text:
using Flux
using Random
using MinimalRLCore

"""
    OnlineAgent

"""
Base.@kwdef mutable struct OnlineAgent{ER<:AbstractReplay, SB, AP<:AbstractValuePolicy, LU, M, TN, O, Φ, SP1, SP2, RT, UC<:Val} <: AbstractAgent
    # Models
    model::M # Currently assumed to be a Flux model
    target_network::TN # Can be either Nothing or a Flux model

    # Learning
    learning_update::LU # Follows the implementation of a Learning Update.
    optimizer::O # Flux optimizer, but could wrap your own (see optimizers.jl)

    # Acting policy (Abstract Value Policy)
    acting_policy::AP # An AbstractValuePolicy. This can be ϵGreedy or something else.

    # State/experience processing.
    replay::ER # Experience Replay. Currently, we only support a vanilla replay but will change in the future based on dispatch.
    state_buffer::SB # A state buffer. Can be Nothing, but the current constructor uses StateBuffer of HistStateBuffer based on hist_length.
    state_preproc::SP1 # A function for processing observations before storage in the state_buffer
    state_postproc::SP2 # Another function for processing the observations after storage in the state_buffer.
    rew_transform::RT # A function for processing the rewards.
    prev_s::Φ # storing the prev_state for adding to the er buffer.

    # params
    batch_size::Int
    target_update_freq::Int
    update_freq::Int
    min_mem_size::Int

    # minor details
    action::Int = 0
    training_steps::Int = 0

end

function OnlineAgent(model,
                     target_network,
                     learning_update,
                     optimizer,
                     acting_policy,
                     replay_size,
                     hist_length,
                     example_state,
                     batch_size,
                     target_update_freq,
                     update_freq,
                     min_mem_size;
                     device = Val{:cpu}(),
                     state_preproc = identity,
                     state_postproc = identity,
                     rew_transform = identity,
                     hist_squeeze=false)

    proc_state = state_preproc(example_state)
    @assert hist_length >= 1
    state_buffer = if hist_length == 1
        DeepRL.StateBuffer{eltype(proc_state)}(replay_size, length(proc_state))
    else
        DeepRL.HistStateBuffer{eltype(proc_state)}(replay_size, size(proc_state), hist_length, hist_squeeze)
    end

    prev_s = if state_buffer isa Nothing
        state_processor(example_state)
    elseif state_buffer isa StateBuffer
        0
    elseif state_buffer isa HistStateBuffer
        zeros(Int, state_buffer.hist_length)
    else
        throw("Unknown StateBuffer please use default Online constructor.")
    end

    replay = ExperienceReplayDef(replay_size, length(prev_s), eltype(prev_s))

    OnlineAgent(model = model,
             target_network = target_network,
             learning_update = learning_update,
             optimizer = optimizer,
             acting_policy = acting_policy,
             replay = replay,
             state_buffer = state_buffer,
             state_preproc = state_preproc,
             state_postproc = state_postproc,
             rew_transform = rew_transform,
             prev_s = prev_s,
             batch_size = batch_size,
             target_update_freq = target_update_freq,
             update_freq = update_freq,
             min_mem_size = min_mem_size,
             device = device)

end

function process_state(agent::OnlineAgent, s; start=false)
    # this stores the state in the state buffer (if not Nothing) and preproccesses the state.
    if agent.state_buffer isa Nothing
        agent.state_preproc(s)
    else
        push!(agent.state_buffer, agent.state_preproc(s); new_episode=start)
        laststate(agent.state_buffer)
    end
end

function get_state_from_buffer(agent::OnlineAgent, s)
    if agent.state_buffer isa Nothing
        to_device(agent.device, agent.state_postproc(s))
    else
        to_device(agent.device, agent.state_postproc(agent.state_buffer[s]))
    end
end

function MinimalRLCore.start!(agent::OnlineAgent,
                              env_s_tp1,
                              rng::AbstractRNG=Random.GLOBAL_RNG)

    # agent.prev_s .= process_state(agent, env_s_tp1)
    agent.prev_s = process_state(agent, env_s_tp1; start=true)
    state = get_state_from_buffer(agent, agent.prev_s)
    resh_state = reshape(state, size(state)..., 1)

    agent.action = sample(agent.acting_policy,
                          to_host(agent.model(resh_state)),
                          rng)

    return agent.action
end

function MinimalRLCore.step!(agent::OnlineAgent,
                             env_s_tp1,
                             r,
                             terminal,
                             rng::AbstractRNG=Random.GLOBAL_RNG)

    proc_state = process_state(agent, env_s_tp1)
    add_ret = add_exp!(agent.replay,
                       (agent.prev_s,
                        findfirst((a)->a==agent.action,
                                  agent.acting_policy.action_set)::Int,
                        proc_state,
                        agent.rew_transform(r), # Atari implementation returns float64s for the reward.
                        terminal))

    update_params!(agent, rng)

    agent.prev_s = proc_state
    prev_s = get_state_from_buffer(agent, agent.prev_s)
    resh_prev_s = reshape(prev_s, size(prev_s)..., 1)

    agent.action = sample(agent.acting_policy,
                          to_host(agent.model(resh_prev_s)),
                          rng)

    return agent.action
end

function update_params!(agent::OnlineAgent, rng)

    if length(agent.replay) > agent.min_mem_size
        if agent.training_steps%agent.update_freq == 0

            e = sample(rng, agent.replay,
                       agent.batch_size)

            s = get_state_from_buffer(agent, e.s)
            r = to_device(agent.device, e.r)
            t = to_device(agent.device, e.t)
            sp = get_state_from_buffer(agent, e.sp)

            ℒ = update!(agent.model,
                        agent.learning_update,
                        agent.optimizer,
                        s,
                        e.a,
                        sp,
                        r,
                        t,
                        agent.target_network)
        end
    end

    # Target network updates
    if !(agent.target_network isa Nothing)
        if agent.training_steps%agent.target_update_freq == 0
            for ps ∈ zip(params(agent.model),
                         params(agent.target_network))
                copyto!(ps[2], ps[1])
            end
        end
    end

    agent.training_steps += 1

    return nothing
end

ext: jl


text:
# Automatically generated using Clang.jl wrap_c, version 0.0.0

using Parameters

@with_kw mutable struct LazGeoKey
    key_id::UInt16 = UInt16(0)
    tiff_tag_location::UInt16 = UInt16(0)
    count::UInt16 = UInt16(0)
    value_offset::UInt16 = UInt16(0)
end

@with_kw mutable struct LazVLR
    reserved::UInt16 = UInt16(0)
    user_id::NTuple{16, UInt8} = ntuple(i -> UInt8(0x0), 16)
    record_id::UInt16 = UInt16(0)
    record_length_after_header::UInt16 = UInt16(0)
    description::NTuple{32, UInt8} = ntuple(i -> UInt8(0x0), 32)
    data::Ptr{UInt8} = pointer("")
end

@with_kw mutable struct LazHeader
    file_source_ID::UInt16 = UInt16(0)
    global_encoding::UInt16 = UInt16(0)
    project_ID_GUID_data_1::UInt32 = UInt32(0)
    project_ID_GUID_data_2::UInt16 = UInt16(0)
    project_ID_GUID_data_3::UInt16 = UInt16(0)
    project_ID_GUID_data_4::NTuple{8, UInt8} = ntuple(i -> UInt8(20), 8)
    # project_ID_GUID_data_4::Array{UInt8, 1} = Array{UInt8, 1}(zeros(0, 8))
    version_major::UInt8 = UInt8(0)
    version_minor::UInt8 = UInt8(0)
    system_identifier::NTuple{32, UInt8} = ntuple(i -> UInt8(20), 32)
    # system_identifier::Array{UInt8} = Array{UInt8}(32)
    generating_software::NTuple{32, UInt8} = ntuple(i -> UInt8(20), 32)
    # generating_software::Array{UInt8, 1} = Array{UInt8, 1}(zeros(0, 32))
    file_creation_day::UInt16 = UInt16(0)
    file_creation_year::UInt16 = UInt16(0)
    header_size::UInt16 = UInt16(0)
    offset_to_point_data::UInt32 = UInt32(0)
    number_of_variable_length_records::UInt32 = UInt32(0)
    point_data_format::UInt8 = UInt8(0)
    point_data_record_length::UInt16 = UInt16(0)
    number_of_point_records::UInt32 = UInt32(0)
    number_of_points_by_return::NTuple{5, UInt32} = ntuple(i -> UInt32(0), 5)
    # number_of_points_by_return::Array{UInt32, 1} = Array{UInt32, 1}(zeros(0, 5))
    x_scale_factor::Float64 = Float64(0.0)
    y_scale_factor::Float64 = Float64(0.0)
    z_scale_factor::Float64 = Float64(0.0)
    x_offset::Float64 = Float64(0.0)
    y_offset::Float64 = Float64(0.0)
    z_offset::Float64 = Float64(0.0)
    max_x::Float64 = Float64(0.0)
    min_x::Float64 = Float64(0.0)
    max_y::Float64 = Float64(0.0)
    min_y::Float64 = Float64(0.0)
    max_z::Float64 = Float64(0.0)
    min_z::Float64 = Float64(0.0)
    start_of_waveform_data_packet_record::UInt64  = UInt64(0)
    start_of_first_extended_variable_length_record::UInt64 = UInt64(0)
    number_of_extended_variable_length_records::UInt32 = UInt32(0)
    extended_number_of_point_records::UInt64 = UInt64(0)
    extended_number_of_points_by_return::NTuple{15, UInt64} = ntuple(i -> UInt64(0), 15)
    # extended_number_of_points_by_return::Array{UInt64, 1} = Array{UInt64, 1}(zeros(0, 15))
    user_data_in_header_size::UInt32 = UInt32(0)
    user_data_in_header::Ptr{UInt8} = pointer("")
    vlrs::Ptr{Cvoid} = pointer("")
    user_data_after_header_size::UInt32 = UInt32(0)
    user_data_after_header::Ptr{UInt8} = pointer("")
end

@with_kw mutable struct LazPoint
    X::Int32 = Int32(0)
    Y::Int32 = Int32(0)
    Z::Int32 = Int32(0)
    intensity::UInt16 = UInt16(0)

    # Evil types
    # U8 return_number : 3;
    # U8 number_of_returns : 3;
    # U8 scan_direction_flag : 1;
    # U8 edge_of_flight_line : 1;  8
    # U8 classification : 5;
    # U8 synthetic_flag : 1;
    # U8 keypoint_flag  : 1;
    # U8 withheld_flag : 1;      8
    return_number::UInt8 = UInt8(0)
    # number_of_returns::UInt8 = UInt8(0)
    # scan_direction_flag::UInt8 = UInt8(0)
    # edge_of_flight_line::UInt8 = UInt8(0)
    classification::UInt8 = UInt8(0)
    # synthetic_flag::UInt8 = UInt8(0)
    # keypoint_flag::UInt8 = UInt8(0)
    # withheld_flag::UInt8 = UInt8(0)
    scan_angle_rank::Int8 = Int8(0)
    user_data::UInt8 = UInt8(0)
    point_source_ID::UInt16 = UInt16(0)

    # Another evil type
    # I16 extended_scan_angle;
    # U8 extended_point_type : 2;
    # U8 extended_scanner_channel : 2;
    # U8 extended_classification_flags : 4;  8
    # U8 extended_classification;
    # U8 extended_return_number : 4;
    # U8 extended_number_of_returns : 4;  8
    extended_scan_angle::Int16 = Int16(0)
    extended_point_type::UInt8 = UInt8(0)
    # extended_scanner_channel::UInt8 = UInt8(0)
    # extended_classification_flags::UInt8 = UInt8(0)
    extended_classification::UInt8 = UInt8(0)
    extended_return_number::UInt8 = UInt8(0)
    # extended_number_of_returns::UInt8 = UInt8(0)
    dummy::NTuple{7, UInt8} = ntuple(i -> UInt8(0), 7)
    gps_time::Float64 = Float64(0.0)
    rgb::NTuple{4, UInt16} = ntuple(i -> UInt16(0), 4)
    wave_packet::NTuple{29, UInt8} = ntuple(i -> UInt8(0), 29)
    num_extra_bytes::Int32 = Int32(0)
    extra_bytes::Ptr{UInt8} = pointer("")
end

return_number(p::LazIO.LazPoint) = (p.return_number & 0b00000111)
number_of_returns(p::LazIO.LazPoint) = (p.return_number & 0b00111000) >> 3
scan_direction(p::LazIO.LazPoint) = Bool((p.return_number & 0b01000000) >> 6)
edge_of_flight_line(p::LazIO.LazPoint) = Bool((p.return_number & 0b10000000) >> 7)


text:
# device management and properties

"""
    device(A::AbstractArray)

Gets the device associated to the Array `A`
"""
device(A::AbstractArray) = error("Not implemented") # COV_EXCL_LINE

"""
Hardware threads of device
"""
threads(device) = error("Not implemented") # COV_EXCL_LINE


text:
using LinearAlgebra: Cholesky

"""
    |(g::GP, c::Observation)

Condition `g` on observation `c`.
"""
function |(g::AbstractGP, c::Observation)
    C = cholesky(Symmetric(cov(c.f)))
    return CompositeGP((|, g, C, C \ (c.y - mean(c.f)), c.f, c.y), g.gpc)
end

const cond_data = Tuple{typeof(|), AbstractGP, Cholesky, AV{<:Real}, FiniteGP, AV{<:Real}}

function mean_vector((_, g, _, α, fx, _)::cond_data, x::AV)
    return mean_vector(g, x) + cov(g, fx.f, x, fx.x) * α
end

function cov((_, g, C_ff, _, fx, _)::cond_data, x::AV)
    return cov(g, x) - Xt_invA_X(C_ff, cov(fx.f, g, fx.x, x))
end
function cov_diag((_, g, C, _, fx, _)::cond_data, x::AV)
    return cov_diag(g, x) - diag_Xt_invA_X(C, cov(fx.f, g, fx.x, x))
end

function cov((_, g, C_ff, _, fx, _)::cond_data, x::AV, x′::AV)
    X, Y = cov(fx.f, g, fx.x, x), cov(fx.f, g, fx.x, x′)
    return cov(g, x, x′) - Xt_invA_Y(X, C_ff, Y)
end
function cov_diag((_, g, C_ff, _, fx, _)::cond_data, x::AV, x′::AV)
    X, Y = cov(fx.f, g, fx.x, x), cov(fx.f, g, fx.x, x′)
    return cov_diag(g, x, x′) - diag_Xt_invA_Y(X, C_ff, Y)
end

function cov((_, g, C, _, fx, _)::cond_data, f′::AbstractGP, x::AV, x′::AV)
    X, Y = cov(fx.f, g, fx.x, x), cov(fx.f, f′, fx.x, x′)
    return cov(g, f′, x, x′) - Xt_invA_Y(X, C, Y)
end

function cov(f::AbstractGP, (_, g, C, _, fx, _)::cond_data, x::AV, x′::AV)
    X, Y = cov(fx.f, f, fx.x, x), cov(fx.f, g, fx.x, x′)
    return cov(f, g, x, x′) - Xt_invA_Y(X, C, Y)
end

function cov_diag((_, g, C, _, fx, _)::cond_data, f′::AbstractGP, x::AV, x′::AV)
    X, Y = cov(fx.f, g, fx.x, x), cov(fx.f, f′, fx.x, x′)
    return cov_diag(g, f′, x, x′) - diag_Xt_invA_Y(X, C, Y)
end

function cov_diag(f::AbstractGP, (_, g, C, _, fx, _)::cond_data, x::AV, x′::AV)
    X, Y = cov(fx.f, f, fx.x, x), cov(fx.f, g, fx.x, x′)
    return cov_diag(f, g, x, x′) - diag_Xt_invA_Y(X, C, Y)
end

function merge(fs::Tuple{Vararg{FiniteGP}})
    block_gp = cross([map(f->f.f, fs)...])
    block_x = BlockData([map(f->f.x, fs)...])
    block_Σy = block_diagonal([map(f->f.Σy, fs)...])
    return FiniteGP(block_gp, block_x, block_Σy)
end
function merge(c::Tuple{Vararg{Observation}})
    block_y = Vector(BlockVector([map(get_y, c)...]))
    return merge(map(get_f, c))←block_y
end
Observation(c::Tuple{Vararg{Observation}}) = merge(c)
Observation(c::Observation...) = Observation(c)

# Multi-arg stuff
|(g::AbstractGP, c::Tuple{Vararg{Observation}}) = g | merge(c)

|(g::Tuple{Vararg{AbstractGP}}, c::Tuple{Vararg{Observation}}) = map(g_ -> g_ | c, g)

|(g::Tuple{Vararg{AbstractGP}}, c::Observation) = g | (c,)


text:
# Note that this script can accept some limited command-line arguments, run
# `julia build_tarballs.jl --help` to see a usage message.
using BinaryBuilder

name     = "CGAL"
rversion = "5.2.1"
version  = VersionNumber(rversion)

# Collection of sources required to build CGAL
sources = [
    ArchiveSource("https://github.com/CGAL/cgal/releases/download/v$rversion/CGAL-$rversion-library.tar.xz",
                  "390aa87c4f21609c19397c4b14abb5ccad3edd1e33933a0089b266f67ce7b111"),
]

# Bash recipe for building across all platforms
script = raw"""
## pre-build setup
# exit on error
set -eu

cmake -B build \
  `# cmake specific` \
  -DCMAKE_BUILD_TYPE=Release \
  -DCMAKE_FIND_ROOT_PATH=$prefix \
  -DCMAKE_INSTALL_PREFIX=$prefix \
  -DCMAKE_TOOLCHAIN_FILE=$CMAKE_TARGET_TOOLCHAIN \
  `# cgal specific` \
  -DCGAL_HEADER_ONLY=OFF \
  -DWITH_CGAL_Core=ON \
  -DWITH_CGAL_ImageIO=ON \
  -DWITH_CGAL_Qt5=OFF \
  CGAL-*/

## and away we go..
cmake --build build --config Release --target install -- -j$nproc
install_license CGAL-*/LICENSE*
"""

# These are the platforms we will build for by default, unless further
# platforms are passed in on the command line
platforms = expand_cxxstring_abis(supported_platforms())

# The products that we will ensure are always built
products = [
    LibraryProduct("libCGAL", :libCGAL),
    LibraryProduct("libCGAL_Core", :libCGAL_Core),
    LibraryProduct("libCGAL_ImageIO", :libCGAL_ImageIO),
]

# Dependencies that must be installed before this package can be built
dependencies = [
    Dependency("boost_jll"),
    Dependency("GMP_jll", v"6.1.2"),
    Dependency("MPFR_jll", v"4.0.2"),
    Dependency("Zlib_jll"),
]

# Build the tarballs.
build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies; preferred_gcc_version = v"7")

max_stars_repo_name: waralex/Yggdrasil

max_issues_repo_name: waralex/Yggdrasil

max_forks_repo_name: waralex/Yggdrasil


text:
"""
    GuidingVectorSampling

Wrapper over any [`AbstractHamiltonian`](@ref) that implements guided vector a.k.a. guided
wave function sampling. In this importance sampling scheme the Hamiltonian is modified as
follows.

```math
\\tilde{H}_{ij} = v_i H_{ij} v_j^{-1}
```

and where `v` is the guiding vector. `v_i` and `v_j` are also thresholded to avoid dividing
by zero (see below).

# Constructors

* `GuidingVectorSampling(::AbstractHamiltonian, vector, eps)`
* `GuidingVectorSampling(::AbstractHamiltonian; vector, eps)`

`eps` is a thresholding parameter used to avoid dividing by zero; all values below `eps` are
set to `eps`. It is recommended that `eps` is in the same value range as the guiding
vector. The default value is set to `eps=norm(v, Inf) * 1e-2`

After construction, we can access the underlying hamiltonian with `G.hamiltonian`, the
`eps` parameter with `G.eps`, and the guiding vector with `G.vector`.

# Example

```jldoctest
julia> H = HubbardReal1D(BoseFS((1,1,1)); u=6.0, t=1.0);

julia> v = DVec(starting_address(H) => 10; capacity=1);

julia> G = GuidingVectorSampling(H, v, 0.1);

julia> get_offdiagonal(H, starting_address(H), 4)
(BoseFS{3,3}((2, 0, 1)), -1.4142135623730951)

julia> get_offdiagonal(G, starting_address(G), 4)
(BoseFS{3,3}((2, 0, 1)), -0.014142135623730952)
```
"""
struct GuidingVectorSampling{A,T,H<:AbstractHamiltonian{T},D,E} <: AbstractHamiltonian{T}
    # The A parameter sets whether this is an adjoint or not.
    # The E parameter is the epsilon value.
    hamiltonian::H
    vector::D
end

function GuidingVectorSampling(h, v::AbstractDVec, eps=1e-2 * norm(v, Inf))
    return GuidingVectorSampling{false,eltype(h),typeof(h),typeof(v),eps}(h, v)
end
function GuidingVectorSampling(h; vector, eps=1e-2 * norm(vector, Inf))
    return GuidingVectorSampling(h, vector, eps)
end

starting_address(h::GuidingVectorSampling) = starting_address(h.hamiltonian)
LOStructure(::Type{<:GuidingVectorSampling{<:Any,<:Any,H}}) where {H} = _lo_str(LOStructure(H))

function LinearAlgebra.adjoint(h::GuidingVectorSampling{A,T,<:Any,D,E}) where {A,T,D,E}
    h_adj = h.hamiltonian'
    return GuidingVectorSampling{!A,T,typeof(h_adj),D,E}(h_adj, h.vector)
end

dimension(::Type{T}, h::GuidingVectorSampling) where T = dimension(T, h.hamiltonian)

Base.getproperty(h::GuidingVectorSampling, s::Symbol) = getproperty(h, Val(s))
Base.getproperty(h::GuidingVectorSampling{<:Any,<:Any,<:Any,<:Any,E}, ::Val{:eps}) where E = E
Base.getproperty(h::GuidingVectorSampling, ::Val{:hamiltonian}) = getfield(h, :hamiltonian)
Base.getproperty(h::GuidingVectorSampling, ::Val{:vector}) = getfield(h, :vector)

# Forward some interface functions.
num_offdiagonals(h::GuidingVectorSampling, add) = num_offdiagonals(h.hamiltonian, add)
diagonal_element(h::GuidingVectorSampling, add) = diagonal_element(h.hamiltonian, add)

_apply_eps(x, eps) = ifelse(iszero(x), eps, ifelse(abs(x) < eps, sign(x) * eps, x))

function guided_vector_modify(value, is_adjoint, eps, guide1, guide2)
    if iszero(guide1) && iszero(guide2)
        return value
    else
        guide1 = _apply_eps(guide1, eps)
        guide2 = _apply_eps(guide2, eps)
        if is_adjoint
            return value * (guide1 / guide2)
        else
            return value * (guide2 / guide1)
        end
    end
end

function get_offdiagonal(h::GuidingVectorSampling{A}, add1, chosen) where A
    add2, matrix_element = get_offdiagonal(h.hamiltonian, add1, chosen)
    guide1 = h.vector[add1]
    guide2 = h.vector[add2]
    return add2, guided_vector_modify(matrix_element, A, h.eps, guide1, guide2)
end

struct GuidingVectorOffdiagonals{
    F,T,A,E,H<:AbstractHamiltonian,D,N<:AbstractOffdiagonals{F,T}
}<:AbstractOffdiagonals{F,T}
    hamiltonian::H
    vector::D
    guide::T
    offdiagonals::N
end

function offdiagonals(h::GuidingVectorSampling{A,T,H,D,E}, a) where {A,T,H,D,E}
    hps = offdiagonals(h.hamiltonian, a)
    guide = h.vector[a]
    return GuidingVectorOffdiagonals{typeof(a),T,A,E,H,D,typeof(hps)}(
        h.hamiltonian, h.vector, guide, hps
    )
end

function Base.getindex(h::GuidingVectorOffdiagonals{F,T,A,E}, i)::Tuple{F,T} where {F,T,A,E}
    add2, matrix_element = h.offdiagonals[i]
    guide1 = h.guide
    guide2 = h.vector[add2]
    return add2, guided_vector_modify(matrix_element, A, E, guide1, guide2)
end

Base.size(h::GuidingVectorOffdiagonals) = size(h.offdiagonals)


text:
type DelayedQL<:FiniteHorizonAgent
    policy::Matrix{UInt64}

    rewards::Array{Float64, 3}

    b::Array{UInt64, 3}
    l::Array{Int64, 3}
    learn::Array{Bool, 3}
    tstar::UInt64
    tval::UInt64

    Q::Array{Float64, 3}
    AU::Array{Float64, 3}

    δ::Float64
    ϵ::Float64
    maxRet::Float64
    mfactor::Float64
    immediate_updates::Bool
end

function DelayedQL(S, A, H, rewards, δ, ϵ, mprecond=1.)
    policy = zeros(UInt64, S, H+1)
    b = zeros(UInt64, S, A, H)
    l = zeros(UInt64, S, A, H)
    AU = zeros(Float64, S, A, H)
    learn = zeros(Bool, S, A, H)
    maxrs = maximum(reshape(rewards, S*A, H), 1)
    V_max = sum(maxrs)
    Q = zeros(Float64, S, A, H+1)
    for t=1:H
        Q[:,:,t] = sum(maxrs[t:H])
    end
    eps1 = ϵ / 3 / H
    mfactor = mprecond * ( 1 + V_max)^2 / 2 / eps1^2 * log(3 * S * A * H / δ *(1 + S * A * H / eps1 ))
    rand!(policy, 1:A) # Initialize with uniformly random policy
    DelayedQL(policy, rewards, b, l, learn, 0, 0, Q, AU, δ, ϵ, V_max, mfactor, true)
end

function observe!(m::DelayedQL, s, a, t, r, sn)
    m.tval += 1
    if m.b[s,a,t] <= m.tstar
        m.learn[s,a,t] = true
    end
    if m.learn[s,a,t]
        if m.l[s,a,t] == 0
            m.b[s,a,t] = m.tval
        end
        m.l[s,a,t] += 1
        m.AU[s,a,t] += m.Q[sn,m.policy[sn,t+1], t+1] + m.rewards[s,a,t]
        if m.immediate_updates
            attempt_update!(m, s,a,t)
        end
    end
end

function attempt_update!(m::DelayedQL, s, a, t)
    if m.l[s,a,t] >= m.mfactor
        if m.Q[s,a,t] >= m.AU[s,a,t] / m.l[s,a,t] + 2 * ϵ1(m)
            m.Q[s,a,t] = m.AU[s,a,t] / m.l[s,a,t] + ϵ1(m)
            m.policy[s,t] = indmax(m.Q[s, :, t])
            m.tstar = m.tval
        elseif m.b[s,a,t] > m.tstar
            m.learn[s,a,t] = false
        end
        m.l[s,a,t] = 0
        m.AU[s,a,t] = 0
    end
end

ϵ1(m::DelayedQL) = m.ϵ / 3 / horizon(m)
update_policy!(m::DelayedQL) = nothing


text:
using Coverage
cov_res = process_folder()
Coveralls.submit(cov_res)
Codecov.submit(cov_res)


text:
using BinaryBuilder, Pkg

name = "AMGX"
version = v"2.1.0"
sources = [
    ArchiveSource("https://github.com/NVIDIA/AMGX/archive/v2.1.0.tar.gz", "6245112b768a1dc3486b2b3c049342e232eb6281a6021fffa8b20c11631f63cc"),
    DirectorySource("./bundled")
]

script = raw"""
# nvcc writes to /tmp, which is a small tmpfs in our sandbox.
# make it use the workspace instead
export TMPDIR=${WORKSPACE}/tmpdir
mkdir ${TMPDIR}

# the build system doesn't find libgcc and libstdc++
if [[ "${nbits}" == 32 ]]; then
    export CFLAGS="-Wl,-rpath-link,/opt/${target}/${target}/lib"
elif [[ "${target}" != *-apple-* ]]; then
    export CFLAGS="-Wl,-rpath-link,/opt/${target}/${target}/lib64"
fi

cd ${WORKSPACE}/srcdir/AMGX*
install_license LICENSE
atomic_patch -p1 "${WORKSPACE}/srcdir/patches/fix-regex-syntax-cmake.patch"

mkdir build
cd build
CMAKE_POLICY_DEFAULT_CMP0021=OLD \
CUDA_BIN_PATH=${prefix}/cuda/bin \
CUDA_LIB_PATH=${prefix}/cuda/lib64 \
CUDA_INC_PATH=${prefix}/cuda/include \
cmake -DCMAKE_TOOLCHAIN_FILE="${CMAKE_TARGET_TOOLCHAIN}" \
      -DCMAKE_INSTALL_PREFIX=${prefix} \
      -DCMAKE_FIND_ROOT_PATH="${prefix}" \
      -DCMAKE_BUILD_TYPE=Release \
      -DCUDA_TOOLKIT_ROOT_DIR="${prefix}/cuda" \
      -DCMAKE_EXE_LINKER_FLAGS="-Wl,-rpath-link,/opt/${target}/${target}/lib64" \
      -Wno-dev \
      ..

make -j${nproc} all
make install

# clean-up
## unneeded static libraries
rm ${libdir}/*.a ${libdir}/sublibs/*.a
"""

platforms = [
    Platform("x86_64", "linux"; libc="glibc", cxxstring_abi = "cxx11"),
]

products = [
    LibraryProduct("libamgxsh", :libamgxsh),
]

dependencies = [
    BuildDependency(PackageSpec(name = "CUDA_full_jll", version = "10.0")),
    Dependency(PackageSpec(name = "CUDA_jll", version = "10.0")),
]

build_tarballs(ARGS, name, version, sources, script, platforms, products, dependencies)

max_stars_repo_name: sharanry/Yggdrasil

max_issues_repo_name: sharanry/Yggdrasil

max_forks_repo_name: sharanry/Yggdrasil


text:
struct ResizableViewVector{T,Data<:AbstractVector{T}} <: AbstractVector{T}
    data::Data
    offset::Int
end

function resizableview(data::AbstractVector, offset::Integer)
    @argcheck firstindex(data) - 1 <= offset <= lastindex(data)
    return ResizableViewVector(data, offset)
end

resizabletail(data::AbstractVector) = resizableview(data, lastindex(data))

Base.size(xs::ResizableViewVector) = (length(xs.data) - xs.offset,)

Base.@propagate_inbounds Base.getindex(xs::ResizableViewVector, i::Int) =
    xs.data[xs.offset+i]

Base.@propagate_inbounds Base.setindex!(xs::ResizableViewVector, x, i::Int) =
    xs.data[xs.offset+i] = x

function Base.resize!(xs::ResizableViewVector, n::Integer)
    resize!(xs.data, xs.offset + n)
    return xs
end

function Base.push!(xs::ResizableViewVector, x)
    push!(xs.data, x)
    return xs
end


text:
using Pkg
using FilePathsBase
using TOML

module Pkgx

export add

function add(pkg::String; bindir=joinpath(home(), ".local", "bin"))
    # Start off adding the package to a temp environment
    dir = mktempdir()
    Pkg.activate(dir)
    Pkg.add(pkg)

    # Check for an "Application.toml" file
    appconfig = joinpath(parent(parent(Path(Base.find_package(pkg)))), "Application.toml")
    if exists(appconfig)
        entrypoints = TOML.parsefile(string(appconfig))["entrypoints"]
    else
        entrypoints = Dict{String,String}()
        entrypoints[lowercase(pkg)] = "julia_main"
    end #if

    # Create a place for the environment to live permanently
    dest = mkpath(joinpath(home(), ".local", "share", "Pkgx", "environments", pkg))
    cp(joinpath(dir, "Project.toml"), dest)
    cp(joinpath(dir, "Manifest.toml"), dest)

    # Create shims
    for pt in entrypoints
        execname = pt.first
        execpath = joinpath(pt.first, execname)
        funname = pt.second
        shimtext =
        """
        #!/bin/sh
        julia --project=$dest -e 'using $pkg; $pkg.$funname()'
        """
        open(string(execpath)) do f
            write(f, shimtext)
        end #do
        chmod(execpath, "+x")
    end #for

end #function

end


text:
module Exp
include("./observable.jl")

experiments = Array{Dict{String,Array{Float64,1}},1}(undef, length(observables))
error_bars = Array{Dict{String,Array{Float64,1}},1}(undef, length(observables))

experiments[observables_index("Nuclear_NFkB")] = Dict(
    "WT" => [
        11.505
        51.52262857
        59.289
        83.68132857
        64.46
        45.027
        34.62007143
        36.81375
        53.31616
        49.223625
        38.17835714
        39.38095714
        44.80163333
        37.23936667
        30.32244
        36.295
        35.46456
        31.04
        34.05333333
        32.175
        29.0548
        30.96566
        35.023
        28.58318571
        22.54135
    ] ./ 83.68132857
)

function get_timepoint(obs_name::String)::Vector{Float64}
    if obs_name == "Nuclear_NFkB"
        [15.0*i for i in 0:24]  # 0, 15, 30, ..., 360 [min.]
    end
end
end # module

max_stars_repo_name: himoto/paramEstim

max_issues_repo_name: himoto/paramEstim


text:
using DiffEqBayes
using Test

@time begin
@time @testset "DynamicHMC" begin include("dynamicHMC.jl") end
@time @testset "Stan_String" begin include("stan_string.jl") end
@time @testset "Stan" begin include("stan.jl") end
@time @testset "Turing" begin include("turing.jl") end # Doesn't work on v0.6
@time @testset "ABC" begin include("abc.jl") end
end


text:
using SlackThreads
using Test
using JET
using Mocking
using JSON3
using CairoMakie
using Logging
using Aqua

Mocking.activate()

const OK_REPLY = Dict("ok" => true)
const FILE_OK_REPLY = Dict("ok" => true, "file" => Dict("permalink" => "LINK"))

function readchomp_reply_patch(reply, count=Ref(0))
    p = @patch function Base.readchomp(cmd)
        count[] += 1
        return JSON3.write(reply)
    end
    return p
end

function readchomp_input_patch(check)
    p = @patch function Base.readchomp(cmd)
        check(cmd)
        return JSON3.write(FILE_OK_REPLY)
    end
    return p
end

function JET_tests(ThreadType; sound_broken=false)
    @testset "JET with args $args" for args in [("hi",), ("hi", "str" => "hello"),
                                                ("hi", "str.txt" => "hello", "a" => "b")]
        thread = withenv(() -> ThreadType(), "SLACK_CHANNEL" => "hi")

        # `@test_call` needs JET v0.5, which needs Julia 1.7
        # We need `@static` here since macroexpansion happens before runtime,
        # i.e. a runtime check is not enough.
        @static if VERSION >= v"1.7"
            @test_call target_modules = (SlackThreads,) thread(args...)
            @test_call target_modules = (SlackThreads,) mode = :sound broken = sound_broken thread(args...)
        end
    end
end

# These can be run in both paths (with exceptions or without),
# since they don't involve throwing exceptions.
function tests_without_errors()
    withenv("SLACK_TOKEN" => "hi", "SLACK_CHANNEL" => nothing) do
        thread = (@test_logs (:warn, r"channel") SlackThread())
        @test thread.ts === nothing
        @test thread.channel === nothing
    end

    withenv("SLACK_TOKEN" => "hi", "SLACK_CHANNEL" => "bye") do
        thread = SlackThread()
        Mocking.apply(readchomp_reply_patch(Dict("ok" => true, "ts" => "abc"))) do
            @test thread("hi").ok == true
            @test thread.ts == "abc"
        end

        @testset "copyto!" begin
            t = SlackThread()
            copyto!(t, (; channel="c", ts="123"))
            @test t.channel == "c"
            @test t.ts == "123"

            # Test we can roundtrip through JSON,
            # then update an existing SlackThread with `copyto!`
            # (since this can be a useful thing to do)
            json = JSON3.write(t)
            new_thread = SlackThread()
            copyto!(new_thread, JSON3.read(json))
            @test new_thread.channel == "c"
            @test new_thread.ts == "123"
        end

        hi_patch = readchomp_input_patch() do cmd
            # Just a reference test; we don't really want to hit up the Slack API
            # from CI here, so let's just check the curl query is one that works
            # from manual testing.
            # This could break for innocuous reasons; in that case, just update it here
            # or find a better test.
            # Currently, this is the only test that checks that the requests we are making
            # are reasonable; we could emit nonsense and all the other tests would pass!
            @test cmd ==
                  `curl -s -X POST -H 'Authorization: Bearer hi' -H 'Content-type: application/json; charset=utf-8' --data '{"channel":"bye","thread_ts":"abc","text":"hi"}' https://slack.com/api/chat.postMessage`
        end

        Mocking.apply(hi_patch) do
            @test thread("hi").ok == true
        end

        option_patch = readchomp_input_patch() do cmd
            # Another reference test
            @test cmd ==
                  `curl -s -X POST -H 'Authorization: Bearer hi' -H 'Content-type: application/json; charset=utf-8' --data '{"channel":"bye","thread_ts":"abc","link_names":true,"text":"hi"}' https://slack.com/api/chat.postMessage`
        end

        Mocking.apply(option_patch) do
            @test thread("hi"; link_names=true).ok == true
        end

        count = Ref(0)
        Mocking.apply(readchomp_reply_patch(Dict("ok" => true, "ts" => "123",
                                                 "file" => Dict("permalink" => "LINK")),
                                            count)) do
            t = SlackThread()
            @test t("bye", "file.txt" => "hi").ok == true
            # `ts` is correctly set
            @test t.ts == "123"
        end
        # we need the general fallback when doing 1 attachment and starting a thread,
        # therefore we get 2 API calls
        @test count[] == 2

        count = Ref(0)
        Mocking.apply(readchomp_reply_patch(OK_REPLY, count)) do
            @test thread("bye").ok == true
            # `ts` doesn't change
            @test thread.ts == "abc"
        end
        @test count[] == 1

        count = Ref(0)
        Mocking.apply(readchomp_reply_patch(FILE_OK_REPLY, count)) do
            @test thread("hi again", "file.txt" => "this is a string",
                         "file2.txt" => "abc").ok == true
        end
        @test count[] == 3 # two file uploads plus the message

        file_patch = readchomp_input_patch() do cmd
            str = string(cmd)
            # Either we're uploading the two files...
            case1 = contains(str, "-F file=@") &&
                    (contains(str, "file.txt") || contains(str, "file2.txt"))
            # Or sending the message with the links
            case2 = contains(str, "hi again<LINK| ><LINK| >")
            @test case1 ⊻ case2
        end

        Mocking.apply(file_patch) do
            @test thread("hi again", "file.txt" => "this is a string",
                         "file2.txt" => "abc").ok == true
        end

        count = Ref(0)
        Mocking.apply(readchomp_reply_patch(OK_REPLY, count)) do
            @test_throws DomainError slack_log_exception(thread) do
                return sqrt(-1)
            end
        end
        @test count[] == 1

        count = Ref(0)
        Mocking.apply(readchomp_reply_patch(OK_REPLY, count)) do
            try
                sqrt(-1)
            catch e
                slack_log_exception(thread, e, catch_backtrace())
            end
        end
        @test count[] == 1

        # Suppress annoying "No strict ticks found" logging messages
        plt1 = with_logger(NullLogger()) do
            return scatter([1.0f0, 2.0f0, 3.0f0], [0.5f0, 4.0f0, 6.0f0])
        end
        plt2 = with_logger(NullLogger()) do
            return scatter(1:10, 1:10)
        end
        Mocking.apply(readchomp_reply_patch(FILE_OK_REPLY)) do
            @test thread("hi", "plot.png" => plt1).ok == true
            @test thread.ts == "abc"

            @test thread("hi", "plot.png" => plt1, "plot2.png" => plt2).ok == true
        end
    end
end

@testset "SlackThreads" begin
    @testset "Aqua" begin
        Aqua.test_all(SlackThreads; ambiguities=false)
    end

    @testset "Utilities" begin
        message_count_suffix = (i, n) -> ""
        messages = SlackThreads.combine_texts(["abcdef", "abc"]; max_length=1,
                                              message_count_suffix)
        @test messages == ["abcdef", "abc"] # 2 messages

        messages = SlackThreads.combine_texts(["abcdef", "abcdef", "a", "b"]; max_length=2,
                                              message_count_suffix)
        @test messages == ["abcdef", "abcdef", "ab"] # 3 messages

        messages = SlackThreads.combine_texts(["abcdef", "abc", "d"]; max_length=4,
                                              message_count_suffix)
        @test messages == ["abcdef", "abcd"]  # can combine last two

        messages = SlackThreads.combine_texts(["d", "abcdef", "abc"]; max_length=4,
                                              message_count_suffix)
        @test messages == ["d", "abcdef", "abc"] # cannot combine anything

        vals = (x for x in ("d", "x", "abcdef", "abc")) # test iterator
        messages = SlackThreads.combine_texts(vals; max_length=4, message_count_suffix)
        @test messages == ["dx", "abcdef", "abc"] # can combine first two

        # Edge case: 1 message
        for max_length in (0, 1, 5)
            messages = SlackThreads.combine_texts(["a"]; max_length, message_count_suffix)
            @test messages == ["a"]
        end

        # Test `message_count_suffix`
        messages = SlackThreads.combine_texts(["a"]; message_count_suffix=(i, n) -> "$i/$n")
        @test messages == ["a1/1"]
    end

    @testset "DummyThreads" begin
        JET_tests(DummyThread)
        d = DummyThread()
        d("hi")
        @test d.logged[1] isa SlackCallRecord
        @test d.logged[1].args == tuple("hi")
        @test d.logged[1].call == :DummyThread
        @test isempty(d.logged[1].kwargs)
        d("hi", "a" => "b")
        @test d.logged[2].args == ("hi", "a" => "b")
        try
            slack_log_exception(d) do
                return sqrt(-1)
            end
        catch DomainError
        end
        @test contains(d.logged[3].args[1], ":alert: Error occured! :alert:")
        @test d.logged[3].call == :send_exception_message

        d("bye"; unfurl_media=true)
        @test d.logged[4].call == :DummyThread
        @test d.logged[4].args == tuple("bye")
        @test d.logged[4].kwargs == (; unfurl_media=true)

        # Properties
        @test propertynames(d) == (:channel, :ts, :logged)
        d.channel = "hi" # no error
        d.ts = "bye" # no error
        @test_throws ErrorException("type DummyThread has no field xyz") d.xyz
        @test_throws MethodError d.channel = 1
        @test_throws MethodError d.ts = 1
        @test d.channel === nothing
        @test d.ts === nothing

        # (de)-Serialization
        roundtrip = JSON3.read(JSON3.write(d), DummyThread)
        @test roundtrip.logged[1] == d.logged[1]
        @test roundtrip.logged[2].kwargs == d.logged[2].kwargs
        # for some reason, JSON3 deserializes a Pair to a Dict (which is a collection of Pairs)
        @test only(roundtrip.logged[2].args[2]) == d.logged[2].args[2]
        @test roundtrip.logged[3] == d.logged[3]
    end

    # Now we test with `SlackThreads.CATCH_EXCEPTIONS[] = false`, i.e.
    # with throwing exceptions. This option exists only for testing, really.
    # The point is we don't want our tests to "pass" while logging exceptions
    # all over the place. If we disable our special "log exceptions instead of throwing"
    # code, do we still pass our tests?
    @testset "With exceptions" begin
        status = SlackThreads.CATCH_EXCEPTIONS[]
        SlackThreads.CATCH_EXCEPTIONS[] = false
        try
            JET_tests(SlackThread; sound_broken=true)

            @testset "Non-throwing tests" begin
                tests_without_errors()
            end

            @testset "`SlackThread` constructor exceptions" begin
                @test_throws MethodError SlackThread(1)
            end

            @testset "`thread` message exceptions" begin
                withenv("SLACK_TOKEN" => "hi", "SLACK_CHANNEL" => "bye") do
                    thread = SlackThread()
                    thread.ts = "abc"

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => true, "ts" => "x"))) do
                        @test thread("bye").ok == true
                        # `ts` still doesn't change
                        @test thread.ts == "abc"

                        @test thread("hi again", "file.txt" => "this is a string").ok ==
                              true

                        # We didn't pass a `file` back in our reply
                        @test_throws KeyError thread("hi again",
                                                     "file.txt" => "this is a string",
                                                     "file2.txt" => "abc")
                    end

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => false,
                                                             "error" => "no"))) do
                        @test_throws SlackThreads.SlackError("no") thread("bye")
                    end

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => false))) do
                        @test_throws SlackThreads.SlackError("No error field returned") thread("bye")
                    end
                end
            end

        finally
            # reset
            SlackThreads.CATCH_EXCEPTIONS[] = status
        end
    end

    # Here, we check that our special "exceptions as logs" machinery works
    # correctly and emits logs.
    @testset "Exceptions as logs" begin
        status = SlackThreads.CATCH_EXCEPTIONS[]
        SlackThreads.CATCH_EXCEPTIONS[] = true
        try
            @testset "Non-throwing tests" begin
                tests_without_errors()
            end

            @testset "`SlackThread` constructor exceptions" begin
                msg = "Error when constructing `SlackThread`"
                thread = (@test_logs (:error, msg) SlackThread(1))
                @test thread.ts === nothing
                @test thread.channel === nothing
            end

            @testset "`thread` message exceptions" begin
                withenv("SLACK_TOKEN" => "hi", "SLACK_CHANNEL" => "bye") do
                    thread = SlackThread()
                    thread.ts = "abc"

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => true, "ts" => "x"))) do
                        @test thread("bye").ok == true
                        # `ts` still doesn't change
                        @test thread.ts == "abc"

                        @test thread("hi again", "file.txt" => "this is a string").ok ==
                              true

                        # We didn't pass a `file` back in our reply
                        msg = "Error when attempting to send message to Slack thread"
                        result = @test_logs (:error, msg) thread("hi again",
                                                                 "file.txt" => "this is a string",
                                                                 "file2.txt" => "abc")
                        @test result === nothing
                    end

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => false,
                                                             "error" => "no"))) do
                        @test (@test_logs (:error, r"Slack API") thread("bye")) === nothing
                    end

                    Mocking.apply(readchomp_reply_patch(Dict("ok" => false))) do
                        @test (@test_logs (:error, r"Slack API") thread("bye")) === nothing
                    end

                    count = Ref(0)
                    Mocking.apply(readchomp_reply_patch(Dict("ok" => false), count)) do
                        @test_logs (:error, "Error reported by Slack API") begin
                            @test_throws DomainError slack_log_exception(thread) do
                                return sqrt(-1)
                            end
                        end
                    end
                    @test count[] == 1
                end
            end
        finally
            # reset
            SlackThreads.CATCH_EXCEPTIONS[] = status
        end
    end
end


text:
import Base: <, iszero, ==, <=, >=, sign, signbit, isinteger, isfinite
export <, iszero, ==, <=, >=, sign, signbit, isinteger, isfinite

function isinteger(a::NumberInterval)
    if floor(a) ≺ a ≺ ceil(a)
        return false
    elseif issingleton(a)
        #NOTE given the previous result, singletons must contain single integers
        return true
    end
    return missing_or_exception(a)
end

function isfinite(a::NumberInterval)
    # per IEEE standard, intervals cannot contain infinities
    if !isempty(a) && !isnan(a)
        return true
    end
    return missing_or_exception(a)
end

function <=(a::NumberInterval, b::NumberInterval)
    if precedes(a, b)
        return true
    elseif precedes(b, a)
        return false
    end
    return missing_or_exception((a, b))
end

function <(a::NumberInterval, b::NumberInterval)
    if strictprecedes(a, b)
        return true
    elseif strictprecedes(b, a)
        return false
    elseif issingleton(a) && issingleton(b) && a ⊆ b
        return false
    end
    return missing_or_exception((a, b))
end

>=(a::NumberInterval, b::NumberInterval) = b <= a

function iszero(a::NumberInterval)
    if !contains_zero(a)
        return false
    end
    if a ⊆ zero(typeof(a))
        return true
    end
    return missing_or_exception(a)
end

function ==(a::NumberInterval, b::NumberInterval)
    if isnan(a) || isnan(b)
        return false
    elseif isdisjoint(a, b)
        return false
    elseif issingleton(a) && issingleton(b) && a ⊆ b
        return true
    end
    return missing_or_exception((a, b))
end

function sign(a::NumberInterval)
    z = zero(typeof(a))
    if a ≺ z
        return -1
    elseif z ≺ a
        return +1
    elseif a ⊆ z
        return  0
    end
    return missing_or_exception(a)
end

function signbit(a::NumberInterval)
    z = zero(typeof(a))
    if a ≺ z
        return true
    elseif precedes(z, a)
        return false
    end
    return missing_or_exception(a)
end


text:
type Config
    splatters::Array{Regex, 1}
    groupers::Array{Regex, 1}
    ignorers::Array{Regex, 1}
end

const defaultconfig = Config(Regex[g"**/_*"], Regex[], Regex[g"**/.*"])

function Config(path::String)
    !isdir(path) && (path = dirname(path))
    # Search for the topmost .raftconfig file
    while !isempty(path)
        f = joinpath(path, ".raftconfig")
        isfile(f) && return parseconfig(f)

        path = dirname(path)
    end

    # Not found; use default configuration
    return defaultconfig
end


text:
import Base: get, getindex, getkey,
            deepcopy, copy, haskey,
            in, display;

export ConstantFunctionDict, CSPDict, CSP, NQueensCSP, SudokuCSP, ZebraCSP,
    get, getkey, getindex, deepcopy, copy, haskey, in,
    assign, unassign, nconflicts, display, infer_assignment,
    MapColoringCSP, backtracking_search, parse_neighbors,
    AC3, first_unassigned_variable, minimum_remaining_values,
    num_legal_values, unordered_domain_values, least_constraining_values,
    no_inference, forward_checking, maintain_arc_consistency,
    min_conflicts, tree_csp_solver, topological_sort,
    support_pruning, solve_zebra;

#Constraint Satisfaction Problems (CSP)

struct ConstantFunctionDict{V}
    value::V

    function ConstantFunctionDict{V}(val::V) where V
        return new(val);
    end
end

ConstantFunctionDict(val) = ConstantFunctionDict{typeof(val)}(val);

copy(cfd::ConstantFunctionDict) = ConstantFunctionDict{typeof(cfd.value)}(cfd.value);

deepcopy(cfd::ConstantFunctionDict) = ConstantFunctionDict{typeof(cfd.value)}(deepcopy(cfd.value));

mutable struct CSPDict
    dict::Union{Nothing, Dict, ConstantFunctionDict}

    function CSPDict(dictionary::Union{Nothing, Dict, ConstantFunctionDict})
        return new(dictionary);
    end
end

function getindex(dict::CSPDict, key)
    if (typeof(dict.dict) <: ConstantFunctionDict)
        return dict.dict.value;
    else
        return getindex(dict.dict, key);
    end
end

function getkey(dict::CSPDict, key, default)
    if (typeof(dict.dict) <: ConstantFunctionDict)
        return dict.dict.value;
    else
        return getkey(dict.dict, key, default);
    end
end

function get(dict::CSPDict, key, default)
    if (typeof(dict.dict) <: ConstantFunctionDict)
        return dict.dict.value;
    else
        return get(dict.dict, key, default);
    end
end

function haskey(dict::CSPDict, key)
    if (typeof(dict.dict) <: ConstantFunctionDict)
        return true;
    else
        return haskey(dict.dict, key);
    end
end

function in(pair::Pair, dict::CSPDict)
    if (typeof(dict.dict) <: ConstantFunctionDict)
        if (getindex(pair, 2) == dict.dict.value)
            return true;
        else
            return false;
        end
    else
        #Call in() function from dict.jl(0.5)/associative.jl(0.6~nightly).
        return in(pair, dict.dict);
    end
end

abstract type AbstractCSP <: AbstractProblem end;

#=

    CSP is a Constraint Satisfaction Problem implementation of AbstractProblem and AbstractCSP.

    This problem contains an unused initial state field to accommodate the requirements

    of some search algorithms.

=#
mutable struct CSP <: AbstractCSP
	vars::AbstractVector
	domains::CSPDict
	neighbors::CSPDict
	constraints::Function
	initial::Tuple
	current_domains::Union{Nothing, Dict}
	nassigns::Int64

    function CSP(vars::AbstractVector, domains::CSPDict, neighbors::CSPDict, constraints::Function;
                initial::Tuple=(), current_domains::Union{Nothing, Dict}=nothing, nassigns::Int64=0)
        return new(vars, domains, neighbors, constraints, initial, current_domains, nassigns)
    end
end

"""
    assign(problem, key, val, assignment)

Overwrite (if an value exists already) assignment[key] with 'val'.
"""
function assign(problem::T, key, val, assignment::Dict) where {T <: AbstractCSP}
    assignment[key] = val;
    problem.nassigns = problem.nassigns + 1;
    nothing;
end

"""
    unassign(problem, key, val, assignment)

Delete the existing (key, val) pair from 'assignment'.
"""
function unassign(problem::T, key, assignment::Dict) where {T <: AbstractCSP}
    if (haskey(assignment, key))
        delete!(assignment, key);
    end
    nothing;
end

function nconflicts(problem::T, key, val, assignment::Dict) where {T <: AbstractCSP}
    return count(
                (function(second_key)
                    return (haskey(assignment, second_key) &&
                        !(problem.constraints(key, val, second_key, assignment[second_key])));
                end),
                problem.neighbors[key]);
end

function display(problem::T, assignment::Dict) where {T <: AbstractCSP}
    println("CSP: ", problem, " with assignment: ", assignment);
    nothing;
end

function actions(problem::T, state::Tuple) where {T <: AbstractCSP}
    if (length(state) == length(problem.vars))
        return [];
    else
        let
            local assignment = Dict(state);
            local var = problem.vars[findfirst((function(e)
                                        return !haskey(assignment, e);
                                    end), problem.vars)];
            return collect((var, val) for val in problem.domains[var]
                            if nconflicts(problem, var, val, assignment) == 0);
        end
    end
end

function get_result(problem::T, state::Tuple, action::Tuple) where {T <: AbstractCSP}
    return (state..., action);
end

function goal_test(problem::T, state::Tuple) where {T <: AbstractCSP}
    let
        local assignment = Dict(state);
        return (length(assignment) == length(problem.vars) &&
                all((function(key)
                            return nconflicts(problem, key, assignment[key], assignment) == 0;
                        end)
                        ,
                        problem.vars));
    end
end

function goal_test(problem::T, state::Dict) where {T <: AbstractCSP}
    let
        local assignment = deepcopy(state);
        return (length(assignment) == length(problem.vars) &&
                all((function(key)
                            return nconflicts(problem, key, assignment[key], assignment) == 0;
                        end),
                        problem.vars));
    end
end

function path_cost(problem::T, cost::Float64, state1::Tuple, action::Tuple, state2::Tuple) where {T <: AbstractCSP}
    return cost + 1;
end

function support_pruning(problem::T) where {T <: AbstractCSP}
    if (problem.current_domains === nothing)
        problem.current_domains = Dict(collect(Pair(key, collect(problem.domains[key])) for key in problem.vars));
    end
    nothing;
end

function suppose(problem::T, key, val) where {T <: AbstractCSP}
    support_pruning(problem);
    local removals::AbstractVector = collect(Pair(key, a) for a in problem.current_domains[key]
                                            if (a != val));
    problem.current_domains[key] = [val];
    return removals;
end

function prune(problem::T, key, value, removals) where {T <: AbstractCSP}
    local not_removed::Bool = true;
    for (i, element) in enumerate(problem.current_domains[key])
        if (element == value)
            deleteat!(problem.current_domains[key], i);
            not_removed = false;
            break;
        end
    end
    if (not_removed)
        error("Could not find ", value, " in ", problem.current_domains[key], " for key '", key, "' to be removed!");
    end
    if (!(typeof(removals) <: Nothing))
        push!(removals, Pair(key, value));
    end
    nothing;
end

function choices(problem::T, key) where {T <: AbstractCSP}
    if (!(problem.current_domains === nothing))
        return problem.current_domains[key];
    else
        return problem.domains[key];
    end
end

function infer_assignment(problem::T) where {T <: AbstractCSP}
    support_pruning(problem);
    return Dict(collect(Pair(key, problem.current_domains[key][1])
                        for key in problem.vars
                            if (1 == length(problem.current_domains[key]))));
end

function restore(problem::T, removals::AbstractVector) where {T <: AbstractCSP}
    for (key, val) in removals
        push!(problem.current_domains[key], val);
    end
    nothing;
end

function conflicted_variables(problem::T, current_assignment::Dict) where {T <: AbstractCSP}
    return collect(var for var in problem.vars
                    if (nconflicts(problem, var, current_assignment[var], current_assignment) > 0));
end

"""
    AC3(problem)

Solve the given problem by applying the arc-consitency algorithm AC-3 (Fig 6.3) to the
given contraint satisfaction problem. Return a boolean indicating whether every arc in
the problem is arc-consistent.
"""
function AC3(problem::T; queue::Union{Nothing, AbstractVector}=nothing, removals::Union{Nothing, AbstractVector}=nothing) where {T <: AbstractCSP}
    if (typeof(queue) <: Nothing)
        queue = collect((X_i, X_k) for X_i in problem.vars for X_k in problem.neighbors[X_i]);
    end
    support_pruning(problem);
    while (length(queue) != 0)
        local X = popfirst!(queue);    #Remove the first item from queue
        local X_i = getindex(X, 1);
        local X_j = getindex(X, 2);
        if (revise(problem, X_i, X_j, removals))
            if (!haskey(problem.current_domains, X_i))
                return false;
            end
            for X_k in problem.neighbors[X_i]
                if (X_k != X_i)
                    push!(queue, (X_k, X_i));
                end
            end
        end
    end
    return true;
end

function revise(problem::T, X_i, X_j, removals::Union{Nothing, AbstractVector}) where {T <: AbstractCSP}
    local revised::Bool = false;
    for x in deepcopy(problem.current_domains[X_i])
        if (all((function(y)
                    return !problem.constraints(X_i, x, X_j, y);
                end),
                problem.current_domains[X_j]))
            prune(problem, X_i, x, removals);
            revised = true;
        end
    end
    return revised;
end

function first_unassigned_variable(problem::T, assignment::Dict) where {T <: AbstractCSP}
    return getindex(problem.vars, findfirst((function(var)
                        return !haskey(assignment, var);
                    end),
                    problem.vars));
end

function minimum_remaining_values(problem::T, assignment::Dict) where {T <: AbstractCSP}
    return argmin_random_tie(collect(v for v in problem.vars if !haskey(assignment, v)),
                            (function(var)
                                return num_legal_values(problem, var, assignment);
                            end));
end

function num_legal_values(problem::T, var, assignment::Dict) where {T <: AbstractCSP}
    if (!(problem.current_domains === nothing))
        return length(problem.current_domains[var]);
    else
        return count((function(val)
                        return (nconflicts(problem, var, val, assignment) == 0);
                    end),
                    problem.domains[var]);
    end
end

function unordered_domain_values(problem::T, var, assignment::Dict) where {T <: AbstractCSP}
    return choices(problem, var);
end

function least_constraining_values(problem::T, var, assignment::Dict) where {T <: AbstractCSP}
    return sort!(deepcopy(choices(problem, var)),
                lt=(function(val)
                       return nconflicts(problem, var, val, assignment);
                    end));
end

function no_inference(problem::T, var, value, assignment::Dict, removals::Union{Nothing, AbstractVector}) where {T <: AbstractCSP}
    return true;
end

function forward_checking(problem::T, var, value, assignment::Dict, removals::Union{Nothing, AbstractVector}) where {T <: AbstractCSP}
    for B in problem.neighbors[var]
        if (!haskey(assignment, B))
            for b in copy(problem.current_domains[B])
                if (!problem.constraints(var, value, B, b))
                    prune(problem, B, b, removals);
                end
            end
            if (length(problem.current_domains[B]) == 0)
                return false;
            end
        end
    end
    return true;
end

function maintain_arc_consistency(problem::T, var, value, assignment::Dict, removals::Union{Nothing, AbstractVector}) where {T <: AbstractCSP}
    return AC3(problem, queue=collect((X, var) for X in problem.neighbors[var]), removals=removals);
end

function backtrack(problem::T, assignment::Dict;
                    select_unassigned_variable::Function=first_unassigned_variable,
                    order_domain_values::Function=unordered_domain_values,
                    inference::Function=no_inference) where {T <: AbstractCSP}
    if (length(assignment) == length(problem.vars))
        return assignment;
    end
    local var = select_unassigned_variable(problem, assignment);
    for value in order_domain_values(problem, var, assignment)
        if (nconflicts(problem, var, value, assignment) == 0)
            assign(problem, var, value, assignment);
            removals = suppose(problem, var, value);
            if (inference(problem, var, value, assignment, removals))
                result = backtrack(problem, assignment,
                                    select_unassigned_variable=select_unassigned_variable,
                                    order_domain_values=order_domain_values,
                                    inference=inference);
                if (!(typeof(result) <: Nothing))
                    return result;
                end
            end
            restore(problem, removals);
        end
    end
    unassign(problem, var, assignment);
    return nothing;
end

"""
    backtracking_search(problem)

Search the given problem by using the backtracking search algorithm (Fig 6.5) and return the solution
if any are found.
"""
function backtracking_search(problem::T;
                            select_unassigned_variable::Function=first_unassigned_variable,
                            order_domain_values::Function=unordered_domain_values,
                            inference::Function=no_inference) where {T <: AbstractCSP}
    local result = backtrack(problem, Dict(),
                            select_unassigned_variable=select_unassigned_variable,
                                    order_domain_values=order_domain_values,
                                    inference=inference);
    if (!(typeof(result) <: Nothing || goal_test(problem, result)))
        error("BacktrackingSearchError: Unexpected result!")
    end
    return result;
end

function parse_neighbors(neighbors::String; vars::AbstractVector=[])
    local new_dict = Dict();
    for var in vars
        new_dict[var] = [];
    end
    local specs::AbstractVector = collect(map(String, split(spec, [':'])) for spec in split(neighbors, [';']));
    for (A, A_n) in specs
        A = String(strip(A));
        if (!haskey(new_dict, A))
            new_dict[A] = [];
        end
        for B in map(String, split(A_n))
            push!(new_dict[A], B);
            if (!haskey(new_dict, B))
                new_dict[B] = [];
            end
            push!(new_dict[B], A);
        end
    end
    return new_dict;
end

function different_values_constraint(A::T1, a::T2, B::T1, b::T2) where {T1, T2}
    return (a != b);
end

function min_conflicts_value(problem::T, var::String, current_assignment::Dict) where {T <: AbstractCSP}
    return argmin_random_tie(problem.domains[var],
                            (function(val)
                                return nconflicts(problem, var, val, current_assignment);
                            end));
end

function min_conflicts_value(problem::T, var::Int64, current_assignment::Dict) where {T <: AbstractCSP}
    return argmin_random_tie(problem.domains[var],
                            (function(val)
                                return nconflicts(problem, var, val, current_assignment);
                            end));
end

"""
    min_conflicts(problem)

Search the given problem by using the min-conflicts algorithm (Fig. 6.8) and return the solution
if any are found.
"""
function min_conflicts(problem::T; max_steps::Int64=100000) where {T <: AbstractCSP}
    local current::Dict = Dict();
    for var in problem.vars
        val = min_conflicts_value(problem, var, current);
        assign(problem, var, val, current);
    end
    for i in 1:max_steps
        local conflicted::AbstractVector = conflicted_variables(problem, current);
        if (length(conflicted) == 0)
            return current;
        end
        local var = rand(RandomDeviceInstance, conflicted);
        local val = min_conflicts_value(problem, var, current);
        assign(problem, var, val, current);
    end
    return nothing;
end

"""
    tree_csp_solver(problem)

Attempt to solve the given problem using the Tree CSP Solver algorithm (Fig. 6.11) and return
the solution if any are found.
"""
function tree_csp_solver(problem::T) where {T <: AbstractCSP}
    local num_of_vars = length(problem.vars);
    local assignment = Dict();
    local root = problem.vars[1];
    X::AbstractVector, parent_dict::Dict = topological_sort(problem, root);
    for X_j in reverse(X)
        if (!make_arc_consistent(problem, parent_dict, X_j))
            return nothing;
        end
    end
    for X_i in X
        if (length(problem.current_domains[X_i]) == 0)
            return nothing;
        end
        assignment[X_i] = problem.current_domains[X_i][1];
    end
    return assignment;
end

"""
    topological_sort(problem, root)

Return a topological sorted AbstractVector and a dictionary of vertices and their parents as keys and values.
The topological depth first search sort first visits the 'root' vertex, traversing the graph by recursively
calling topological_sort_postorder_dfs on vertices returned by problem.neighbors[vertex].

The dictionary returned does not have keys (vertices) that do not have a parent vertex.
"""
function topological_sort(problem::T, root::String) where {T <: AbstractCSP}
    local sorted_nodes = [];
    local parents = Dict();
    local visited = Dict();
    for v in problem.vars
        visited[v] = false;
    end
    local neighbors = problem.neighbors;
    topological_sort_postorder_dfs(root, visited, neighbors, parents, sorted_nodes, nothing);
    return sorted_nodes, parents;
end

function topological_sort_postorder_dfs(vertex::String,
                                        visited::Dict,
                                        neighbors::CSPDict,
                                        parents::Dict,
                                        sorted_vertices::AbstractVector,
                                        parent::Union{Nothing, String})
    visited[vertex] = true;
    for w in neighbors[vertex]
        if (!visited[w])
            topological_sort_postorder_dfs(w, visited, neighbors, parents, sorted_vertices, vertex);
        end
    end
    insert!(sorted_vertices, 1, vertex);
    if (!(typeof(parent) <: Nothing))
        parents[vertex] = parent;
    end
    nothing;
end

function make_arc_consistent(problem::T, parent::Dict, X_j) where {T <: AbstractCSP}
    println("make_arc_consistent() is not implemented yet for ", typeof(problem), "!");
    nothing;
end

function MapColoringCSP(colors::AbstractVector, neighbors::String)
    local parsed_neighbors = parse_neighbors(neighbors);
    return CSP(collect(keys(parsed_neighbors)), CSPDict(ConstantFunctionDict(colors)), CSPDict(parsed_neighbors), different_values_constraint);
end

function MapColoringCSP(colors::AbstractVector, neighbors::Dict)
    return CSP(collect(keys(neighbors)), CSPDict(ConstantFunctionDict(colors)), CSPDict(neighbors), different_values_constraint);
end

australia_csp = MapColoringCSP(["R", "G", "B"], "SA: WA NT Q NSW V; NT: WA Q; NSW: Q V; T: ");

usa_csp = MapColoringCSP(["R", "G", "B", "Y"],
                        "WA: OR ID; OR: ID NV CA; CA: NV AZ; NV: ID UT AZ; ID: MT WY UT;
                        UT: WY CO AZ; MT: ND SD WY; WY: SD NE CO; CO: NE KA OK NM; NM: OK TX;
                        ND: MN SD; SD: MN IA NE; NE: IA MO KA; KA: MO OK; OK: MO AR TX;
                        TX: AR LA; MN: WI IA; IA: WI IL MO; MO: IL KY TN AR; AR: MS TN LA;
                        LA: MS; WI: MI IL; IL: IN KY; IN: OH KY; MS: TN AL; AL: TN GA FL;
                        MI: OH IN; OH: PA WV KY; KY: WV VA TN; TN: VA NC GA; GA: NC SC FL;
                        PA: NY NJ DE MD WV; WV: MD VA; VA: MD DC NC; NC: SC; NY: VT MA CT NJ;
                        NJ: DE; DE: MD; MD: DC; VT: NH MA; MA: NH RI CT; CT: RI; ME: NH;
                        HI: ; AK: ");

france_csp = MapColoringCSP(["R", "G", "B", "Y"],
                            "AL: LO FC; AQ: MP LI PC; AU: LI CE BO RA LR MP; BO: CE IF CA FC RA
                            AU; BR: NB PL; CA: IF PI LO FC BO; CE: PL NB NH IF BO AU LI PC; FC: BO
                            CA LO AL RA; IF: NH PI CA BO CE; LI: PC CE AU MP AQ; LO: CA AL FC; LR:
                            MP AU RA PA; MP: AQ LI AU LR; NB: NH CE PL BR; NH: PI IF CE NB; NO:
                            PI; PA: LR RA; PC: PL CE LI AQ; PI: NH NO CA IF; PL: BR NB CE PC; RA:
                            AU BO FC PA LR");

function queen_constraint(A, a, B, b)
    return ((A == B) || ((a != b) && (A + a != B + b) && (A - a != B - b)));
end

#=

    NQueensCSP is a N-Queens Constraint Satisfaction Problem implementation of AbstractProblem and AbstractCSP.

=#
mutable struct NQueensCSP <: AbstractCSP
    vars::AbstractVector
    domains::CSPDict
    neighbors::CSPDict
    constraints::Function
    initial::Tuple
    current_domains::Union{Nothing, Dict}
    nassigns::Int64
    rows::AbstractVector
    backslash_diagonals::AbstractVector
    slash_diagonals::AbstractVector

    function NQueensCSP(n::Int64; initial::Tuple=(), current_domains::Union{Nothing, Dict}=nothing, nassigns::Int64=0)
        return new(collect(1:n),
                    CSPDict(ConstantFunctionDict(collect(1:n))),
                    CSPDict(ConstantFunctionDict(collect(1:n))),
                    queen_constraint,
                    initial,
                    current_domains,
                    nassigns,
                    fill(0, n),
                    fill(0, ((2 * n) - 1)),
                    fill(0, ((2 * n) - 1)),);
    end
end

function nconflicts(problem::NQueensCSP, key::Int64, val::Int64, assignment::Dict)
    local num_of_vars::Int64 = length(problem.vars);
    local c::Int64 = problem.rows[val] +
                    problem.backslash_diagonals[key + val - 1] +
                    problem.slash_diagonals[key - val + num_of_vars];
    if (get(assignment, key, nothing) == val)
        c = c - 3;
    end
    return c;
end

function assign(problem::NQueensCSP, key::Int64, val::Int64, assignment::Dict)
    local old_val = get(assignment, key, nothing);
    if (old_val != val)
        if (!(typeof(old_val) <: Nothing))
            record_conflict(problem, assignment, key, old_val, -1);
        end
        record_conflict(problem, assignment, key, val, 1);

        assignment[key] = val;
        problem.nassigns = problem.nassigns + 1;
    end
    nothing;
end

function unassign(problem::NQueensCSP, key::Int64, assignment::Dict)
    if (haskey(assignment, key))
        record_conflict(problem, assignment, key, assignment[key], -1);
        delete!(assignment, key);
    end
    nothing;
end

function record_conflict(problem::NQueensCSP, assignment::Dict, key::Int64, val::Int64, delta::Int64)
    local num_of_vars::Int64 = length(problem.vars);
    problem.rows[val] = problem.rows[val] + delta;
    problem.backslash_diagonals[key + val - 1] = problem.backslash_diagonals[key + val - 1] + delta;
    problem.slash_diagonals[key - val + num_of_vars] = problem.slash_diagonals[key - val + num_of_vars] + delta;
end

function display(problem::NQueensCSP, assignment::Dict)
    local num_of_vars::Int64 = length(problem.vars);
    for val in 1:num_of_vars
        for key in 1:num_of_vars
            local piece::String;
            if (get(assignment, key, "") == val)
                piece = "Q";
            elseif ((key + val - 1) % 2 == 0)
                piece = ".";
            else
                piece = "-";
            end
            print(piece);
        end
        print("    ");
        for key in 1:num_of_vars
            local piece::String;
            if (get(assignment, key, "") == val)
                piece = "*";
            else
                piece = " ";
            end
            print(nconflicts(problem, key, val, assignment), piece);
        end
        println();
    end
    nothing;
end

easy_sudoku_grid = "..3.2.6..9..3.5..1..18.64....81.29..7.......8..67.82....26.95..8..2.3..9..5.1.3..";

harder_sudoku_grid = "4173698.5.3..........7......2.....6.....8.4......1.......6.3.7.5..2.....1.4......";

struct SudokuInitialState
    index_grid::AbstractVector
    boxes::AbstractVector
    rows::AbstractVector
    cols::AbstractVector

    function SudokuInitialState()
        local index_iter = Base.Iterators.countfrom();
        local index_position::AbstractVector = [0];
        local index_grid = collect(collect(collect(collect((function(it, ip)
                                                                ip[1]=iterate(it, ip[1])[2];
                                                                return ip[1];
                                                            end)(index_iter, index_position) for x in 1:3)
                                                    for y in 1:3)
                                            for box_x in 1:3)
                                    for box_y in 1:3);
        local boxes = reduce(vcat, collect(collect(reduce(vcat, row) for row in box_row) for box_row in index_grid));
        local rows = reduce(vcat, collect(collect(reduce(vcat, collect(uneval))
                                                for uneval in zip(box_row...,))
                                        for box_row in index_grid));
        local cols = collect(collect(eval_tuple) for eval_tuple in zip(rows...,));
        return new(index_grid, boxes, rows, cols);
    end
end

sudoku_indices = SudokuInitialState();

#=

    SudokuCSP is a Sudoku Constraint Satisfaction Problem implementation of AbstractProblem and AbstractCSP.

=#
mutable struct SudokuCSP <: AbstractCSP
    vars::AbstractVector
    domains::CSPDict
    neighbors::CSPDict
    constraints::Function
    initial::Tuple
    current_domains::Union{Nothing, Dict}
    nassigns::Int64

    function SudokuCSP(grid::String; initial::Tuple=(), current_domains::Union{Nothing, Dict}=nothing, nassigns::Int64=0)
        local neighbors = Dict{Int64, Any}(collect(Pair(cell, Set()) for cell in reduce(vcat, sudoku_indices.rows)));
        for unit in map(Set, vcat(sudoku_indices.boxes, sudoku_indices.rows, sudoku_indices.cols))
            for cell in unit
                neighbors[cell] = union(neighbors[cell], setdiff(unit, Set(cell)));
            end
        end

        local squares = map(String, collect(m.match for m in eachmatch(r"\d|\.", grid)));
        if (length(squares) != length(reduce(vcat, sudoku_indices.rows)))
            error("SudokuCSPError: Invalid Sudoku grid!")
        end
        local domains = Dict(collect(Pair(key,
                                        if_((number_str in ["1", "2", "3", "4", "5", "6", "7", "8", "9"]),
                                            [Int64(Vector{UInt8}(number_str)[1]) - 48],
                                            [1, 2, 3, 4, 5, 6, 7, 8, 9]))
                                    for (key, number_str) in zip(reduce(vcat, sudoku_indices.rows), squares)));
        #Sort the keys of 'domains' when creating the 'vars' field of the new SudokuCSP problem.
        #Otherwise, backtracking search on SudokuCSP will run indefinitely.
        return new(sort(collect(keys(domains))), CSPDict(domains), CSPDict(neighbors), different_values_constraint,
                    initial, current_domains, nassigns);
    end
end

function display(problem::SudokuCSP, assignment::Dict)
    return join(collect(
                join(reduce(
                            (function(lines1, lines2)
                                return map(
                                            (function(array_str)
                                                return join(array_str, " | ");
                                            end),
                                            zip(lines1, lines2));
                            end),
                            map((function(box)
                                    return collect(join(map(
                                                            (function(cell)
                                                                return string(get(assignment, cell, "."))
                                                            end),
                                                            row), " ")
                                                    for row in box);
                                end),
                                box_row)), "\n")
                    for box_row in sudoku_indices.index_grid),
            "\n------+-------+------\n");
end

struct ZebraInitialState
    colors::Array{String, 1}
    pets::Array{String, 1}
    drinks::Array{String, 1}
    countries::Array{String, 1}
    smokes::Array{String, 1}

    function ZebraInitialState()
        local colors = map(String, split("Red Yellow Blue Green Ivory"));
        local pets = map(String, split("Dog Fox Snails Horse Zebra"));
        local drinks = map(String, split("OJ Tea Coffee Milk Water"));
        local countries = map(String, split("Englishman Spaniard Norwegian Ukranian Japanese"));
        local smokes = map(String, split("Kools Chesterfields Winston LuckyStrike Parliaments"));
        return new(colors, pets, drinks, countries, smokes);
    end
end

zebra_constants = ZebraInitialState();

function zebra_constraint(A::String, a, B::String, b; recursed::Bool=false)
    local same::Bool = (a == b);
    local next_to::Bool = (abs(a - b) == 1);
    if (A == "Englishman" && B == "Red")
        return same;
    elseif (A == "Spaniard" && B == "Dog")
        return same;
    elseif (A == "Chesterfields" && B == "Fox")
        return next_to;
    elseif (A == "Norwegian" && B == "Blue")
        return next_to;
    elseif (A == "Kools" && B == "Yellow")
        return same;
    elseif (A == "Winston" && B == "Snails")
        return same;
    elseif (A == "LuckyStrike" && B == "OJ")
        return same;
    elseif (A == "Ukranian" && B == "Tea")
        return same;
    elseif (A == "Japanese" && B == "Parliaments")
        return same;
    elseif (A == "Kools" && B == "Horse")
        return next_to;
    elseif (A == "Coffee" && B == "Green")
        return same;
    elseif (A == "Green" && B == "Ivory")
        return ((a - 1) == b);
    elseif (!recursed)
        return zebra_constraint(B, b, A, a, recursed=true);
    elseif ((A in zebra_constants.colors && B in zebra_constants.colors) ||
            (A in zebra_constants.pets && B in zebra_constants.pets) ||
            (A in zebra_constants.drinks && B in zebra_constants.drinks) ||
            (A in zebra_constants.countries && B in zebra_constants.countries) ||
            (A in zebra_constants.smokes && B in zebra_constants.smokes))
        return !same;
    else
        error("ZebraConstraintError: This constraint could not be evaluated on the given arguments!");
    end
end

#=

    ZebraCSP is a Zebra Constraint Satisfaction Problem implementation of AbstractProblem and AbstractCSP.

=#
mutable struct ZebraCSP <: AbstractCSP
    vars::AbstractVector
    domains::CSPDict
    neighbors::CSPDict
    constraints::Function
    initial::Tuple
    current_domains::Union{Nothing, Dict}
    nassigns::Int64

    function ZebraCSP(;initial::Tuple=(), current_domains::Union{Nothing, Dict}=nothing, nassigns::Int64=0)
        local vars = vcat(zebra_constants.colors,
                        zebra_constants.pets,
                        zebra_constants.drinks,
                        zebra_constants.countries,
                        zebra_constants.smokes);
        local domains = Dict();
        for var in vars
            domains[var] = collect(1:5);
        end
        domains["Norwegian"] = [1];
        domains["Milk"] = [3];
        neighbors = parse_neighbors("Englishman: Red;
                Spaniard: Dog; Kools: Yellow; Chesterfields: Fox;
                Norwegian: Blue; Winston: Snails; LuckyStrike: OJ;
                Ukranian: Tea; Japanese: Parliaments; Kools: Horse;
                Coffee: Green; Green: Ivory", vars=vars);
        for category in [zebra_constants.colors,
                        zebra_constants.pets,
                        zebra_constants.drinks,
                        zebra_constants.countries,
                        zebra_constants.smokes]
            for A in category
                for B in category
                    if (A != B)
                        if (!(B in neighbors[A]))
                            push!(neighbors[A], B);
                        end
                        if (!(A in neighbors[B]))
                            push!(neighbors[B], A);
                        end
                    end
                end
            end
        end
        return new(vars, CSPDict(domains), CSPDict(neighbors), zebra_constraint, initial, current_domains, nassigns);
    end
end

function solve_zebra(problem::ZebraCSP, algorithm::Function; kwargs...)
    local answer = algorithm(problem; kwargs...);
    for house in collect(1:5)
        print("House ", house);
        for (key, val) in collect(answer)
            if (val == house)
                print(" ", key);
            end
        end
        println();
    end
    return answer["Zebra"], answer["Water"], problem.nassigns, answer;
end

max_stars_repo_name: mikhail-j/aima-julia-aimacode-fork

max_issues_repo_name: mikhail-j/aima-julia-aimacode-fork

max_forks_repo_name: mikhail-j/aima-julia-aimacode-fork