text:
module FastArrayOps
import Base.LinAlg: BlasReal, BlasComplex, BlasFloat, BlasInt, BlasChar
const libblas = Base.libblas_name
export fast_scale!,     unsafe_fast_scale!,
       fast_add!,       unsafe_fast_add!,
       fast_addscal!,   unsafe_fast_addscal!,
       fast_copy!,      unsafe_fast_copy!,
       fast_fill!,      unsafe_fast_fill!
export fast_check1, fast_check2, fast_check3, nmax2nel, nel2nmax, fast_args2range, fast_range2args

# WARNING: FastArrayOps.jl gets overwritten by FastArrayOps_src.jl when running make.jl

## CONSTANTS

const NLIM_SCALE = 13
const NLIM_SCALE_OOP1 = 80
const NLIM_SCALE_OOP2 = 100000
const NLIM_SCALEARR = typemax(Int)
const NLIM_SCALEARR_OOP = typemax(Int)
const NLIM_ADD = typemax(Int)
const NLIM_ADD_OOP = typemax(Int)
const NLIM_ADDARR = 13
const NLIM_ADDARR_OOP1 = 30
const NLIM_ADDARR_OOP2 = 100000
const NLIM_ADDARRSCAL = 13
const NLIM_ADDARRSCAL_OOP1 = 30
const NLIM_ADDARRSCAL_OOP2 = 100000
const NLIM_COPY1 = 80
const NLIM_COPY2 = 100000
const NLIM_FILL = 13

## UTILS

function nmax2nel(i::Int, inc::Int, nmax::Int)
    @assert 0 < i
    nmax < i && return 0
    return div(nmax - i, abs(inc)) + 1
end
function nel2nmax(i::Int, inc::Int, nel::Int)
    @assert 0 < i
    nel < 0 && return i - 1
    return i + (nel- 1)*abs(inc)
end
function fast_args2range(i::Int, inc::Int, n::Int)
    @assert 0 < i
    r = i:abs(inc):nel2nmax(i, abs(inc), n)
    if inc < 0
        r = reverse(r)
    end
    return r
end
function fast_range2args(r::Range)
    inc = step(r)
    if inc > 0
        return (first(r), inc, length(r))
    else
        return (last(r), inc, length(r))
    end
end

function fast_check1(x, ix, incx, n)
    # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    return 0
end
function fast_check2(x, ix, incx, y, iy, incy, n)
    # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    return 0
end
function fast_check3(x, ix, incx, y, iy, incy, z, iz, incz, n)
    # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    return 0
end

# utils for 0 fill
inttype(::Type{Float32}) = Int32
inttype(::Type{Float64}) = Int64
inttype(::Type{Complex64}) = Int64
inttype(::Type{Complex128}) = Int128
function fast_reinterpret1{T<:Number}(::Type{T}, a::Array)
    @assert length(a) == 1
    ccall(:jl_reshape_array, Array{T,1}, (Any, Any, Any), Array{T,1}, a, (1,))
end
# are float 0s zero bits?
function zerobits()
    v = reinterpret(inttype(Float32), convert(Float32, 0))
    v += reinterpret(inttype(Float64), convert(Float64, 0))
    v += fast_reinterpret1(inttype(Complex64), [convert(Complex64, 0)])[1]
    v += fast_reinterpret1(inttype(Complex128), [convert(Complex128, 0)])[1]
    if v == 0
        return true
    else
        return false
    end
end
const ZEROFLOAT = zerobits()

for (fscal, fcopy, faxpy, ftbmv, fsbmv, elty) in (
                        (:dscal_, :dcopy_, :daxpy_, :dtbmv_, :dsbmv_, :Float64),
                        (:sscal_, :scopy_, :saxpy_, :stbmv_, :ssbmv_, :Float32),
                        (:zscal_, :zcopy_, :zaxpy_, :ztbmv_, :zsbmv_, :Complex128),
                        (:cscal_, :ccopy_, :caxpy_, :ctbmv_, :csbmv_, :Complex64))

## SCALE METHODS

for (f, isunsafe) in ( (:fast_scale!, false), (:unsafe_fast_scale!, true) )
@eval begin

# x = a*x
# =======
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    if n < $NLIM_SCALE*incx
        # src: scalarr1_for.jl
incx = abs(incx)
@inbounds for i = ix:incx:ix-1+n*incx
    x[i] = *( a, x[i] )
end
    else
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set1_inc1.jl
incx = 1
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    if n < $NLIM_SCALE #*incx
        # src: scalarr1_for_inc1.jl
@inbounds for i = ix:ix-1+n
    x[i] = *( a, x[i] )
end
    else
        # src: set1_inc1.jl
incx = 1
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end

# x = a*y
# =======
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    mul = max(abs(incx), abs(incy))
    if n < $NLIM_SCALE_OOP1*mul || n*mul > $NLIM_SCALE_OOP2
        # src: scalarr1_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = *( a, y[iy+i*incy] )
end
    else
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    mul = abs(incx)
    if n < $NLIM_SCALE_OOP1*mul || n*mul > $NLIM_SCALE_OOP2
        # src: scalarr1_foroop_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = *( a, y[d+i] )
end
    else
        # src: set2_inceq.jl
incy = incx
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_SCALE_OOP1 || n > $NLIM_SCALE_OOP2
        # src: scalarr1_foroop_inc1.jl
d = iy - ix
@inbounds for i = ix:ix-1+n
    x[i] = *( a, y[d+i] )
end
    else
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_SCALE_OOP1 || n > $NLIM_SCALE_OOP2
        # src: scalarr1_foroop_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = *( a, y[i] )
end
    else
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        # src: blas_scale.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
ccall(($(string(fscal)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), px, &(incx))
    end
    return x
end

# x = x.*y
# ========
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: arr2xy_for.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = *( x[ix+i*incx], y[iy+i*incy] )
end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: arr2xy_for_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = *( x[i], y[d+i] )
end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: arr2xy_for_inc1.jl
d = iy - ix
@inbounds for i = ix:ix-1+n
    x[i] = *( x[i], y[d+i] )
end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: arr2xy_for_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = *( x[i], y[i] )
end
    return x
end

# x = y.*z
# ========

# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, z::Array{$elty}, iz::Int, incz::Int, n::Int)
    $isunsafe || begin
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    # src: arr2yz_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
incz < 0 && (iz = iz+(n-1)*abs(incz))
@inbounds for i = 0:n-1
    x[ix+i*incx] = *( y[iy+i*incy], z[iz+i*incz] )
end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, n::Int)
    $isunsafe || begin
        # src: set3_inceq.jl
incy = incx
incz = incx
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    # src: arr2yz_foroop_inceq.jl
incx = abs(incx)
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = *( y[dy+i], z[dz+i] )
end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, n::Int)
    $isunsafe || begin
        # src: set3_inc1.jl
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    # src: arr2yz_foroop_inc1.jl
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:ix-1+n
    x[i] = *( y[dy+i], z[dz+i] )
end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, z::Array{$elty}, n::Int)
    $isunsafe || begin
        # src: set3_inc1ieq.jl
iy = ix
iz = ix
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    # src: arr2yz_foroop_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = *( y[i], z[i] )
end
    return x
end

end # eval begin
end # for

## ADD METHODS

for (f, isunsafe) in ( (:fast_add!, false), (:unsafe_fast_add!, true) )
@eval begin

# x = x + a
# =======
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    # src: scalarr1_for.jl
incx = abs(incx)
@inbounds for i = ix:incx:ix-1+n*incx
    x[i] = +( a, x[i] )
end
    return x
end

# inc1
function ($f)(x::Array{$elty}, ix::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set1_inc1.jl
incx = 1
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    # src: scalarr1_for_inc1.jl
@inbounds for i = ix:ix-1+n
    x[i] = +( a, x[i] )
end
    return x
end

# x = y + a
# =======
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: scalarr1_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = +( a, y[iy+i*incy] )
end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: scalarr1_foroop_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = +( a, y[d+i] )
end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: scalarr1_foroop_inc1.jl
d = iy - ix
@inbounds for i = ix:ix-1+n
    x[i] = +( a, y[d+i] )
end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: scalarr1_foroop_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = +( a, y[i] )
end
    return x
end

# x = x + y
# =========
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR #*mul # || n*mul > $NLIM_SCALEARR
        # src: arr2xy_for.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = +( x[ix+i*incx], y[iy+i*incy] )
end
    else
        a = 1
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR
        # src: arr2xy_for_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = +( x[i], y[d+i] )
end
    else
        a = 1
        # src: set2_inceq.jl
incy = incx
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR
        # src: arr2xy_for_inc1.jl
d = iy - ix
@inbounds for i = ix:ix-1+n
    x[i] = +( x[i], y[d+i] )
end
    else
        a = 1
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR
        # src: arr2xy_for_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = +( x[i], y[i] )
end
    else
        a = 1
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end

# x = y + z
# =========
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, z::Array{$elty}, iz::Int, incz::Int, n::Int)
    $isunsafe || begin
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR_OOP1 || n > $NLIM_ADDARR_OOP2 #*mul # || n*mul > $NLIM_SCALEARR
        # src: arr2yz_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
incz < 0 && (iz = iz+(n-1)*abs(incz))
@inbounds for i = 0:n-1
    x[ix+i*incx] = +( y[iy+i*incy], z[iz+i*incz] )
end
    else
        a = 1
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, n::Int)
    $isunsafe || begin
        # src: set3_inceq.jl
incy = incx
incz = incx
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR_OOP1 || n > $NLIM_ADDARR_OOP2
        # src: arr2yz_foroop_inceq.jl
incx = abs(incx)
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = +( y[dy+i], z[dz+i] )
end
    else
        a = 1
        # src: set3_inceq.jl
incy = incx
incz = incx
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, n::Int)
    $isunsafe || begin
        # src: set3_inc1.jl
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR_OOP1 || n > $NLIM_ADDARR_OOP2
        # src: arr2yz_foroop_inc1.jl
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:ix-1+n
    x[i] = +( y[dy+i], z[dz+i] )
end
    else
        a = 1
        # src: set3_inc1.jl
incx = 1
incy = 1
incz = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, z::Array{$elty}, n::Int)
    $isunsafe || begin
        # src: set3_inc1ieq.jl
iy = ix
iz = ix
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARR_OOP1 || n > $NLIM_ADDARR_OOP2
        # src: arr2yz_foroop_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = +( y[i], z[i] )
end
    else
        a = 1
        # src: set3_inc1ieq.jl
iy = ix
iz = ix
incx = 1
incy = 1
incz = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end

end # eval begin
end # for

## ADDSCAL METHODS

for (f, isunsafe) in ( (:fast_addscal!, false), (:unsafe_fast_addscal!, true) )
@eval begin

# x = x + a*y
# =========
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL #*mul # || n*mul > $NLIM_SCALEARR
        # src: addarrscal_for.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = x[ix+i*incx] + y[iy+i*incy]*a
end
    else
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL
        # src: addarrscal_for_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = x[i] + y[d+i]*a
end
    else
        # src: set2_inceq.jl
incy = incx
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL
        # src: addarrscal_for_inc1.jl
d = iy - ix
@inbounds for i = ix:ix-1+n
    x[i] = x[i] + y[d+i]*a
end
    else
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, a::$elty, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL
        # src: addarrscal_for_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = x[i] + y[i]*a
end
    else
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end

# x = y + a*z
# =========
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, z::Array{$elty}, iz::Int, incz::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL_OOP1 || n > $NLIM_ADDARRSCAL_OOP2 #*mul # || n*mul > $NLIM_SCALEARR
        # src: addarrscal_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
incz < 0 && (iz = iz+(n-1)*abs(incz))
@inbounds for i = 0:n-1
    x[ix+i*incx] = y[iy+i*incy] + z[iz+i*incz]*a
end
    else
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set3_inceq.jl
incy = incx
incz = incx
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL_OOP1 || n > $NLIM_ADDARRSCAL_OOP2
        # src: addarrscal_foroop_inceq.jl
incx = abs(incx)
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = y[dy+i] + z[dz+i]*a
end
    else
        # src: set3_inceq.jl
incy = incx
incz = incx
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, z::Array{$elty}, iz::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set3_inc1.jl
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL_OOP1 || n > $NLIM_ADDARRSCAL_OOP2
        # src: addarrscal_foroop_inc1.jl
dy = iy - ix
dz = iz - ix
@inbounds for i = ix:ix-1+n
    x[i] = y[dy+i] + z[dz+i]*a
end
    else
        # src: set3_inc1.jl
incx = 1
incy = 1
incz = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, z::Array{$elty}, a::$elty, n::Int)
    $isunsafe || begin
        # src: set3_inc1ieq.jl
iy = ix
iz = ix
incx = 1
incy = 1
incz = 1
        # src: fast_check3.jl
(0 != incx && 0 != incy && 0 != incz) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy && 0 < iz) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
iz+(n-1)*abs(incz) <= length(z) || throw(BoundsError())
    end
    if n < $NLIM_ADDARRSCAL_OOP1 || n > $NLIM_ADDARRSCAL_OOP2
        # src: addarrscal_foroop_inc1ieq.jl
@inbounds for i = ix:ix-1+n
    x[i] = y[i] + z[i]*a
end
    else
        # src: set3_inc1ieq.jl
iy = ix
iz = ix
incx = 1
incy = 1
incz = 1
        # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
        y, iy, incy = z, iz, incz
        # src: blas_axpy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(faxpy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), &(a), py, &(incy), px, &(incx))
    end
    return x
end

end # eval begin
end # for

## COPY METHODS

for (f, isunsafe) in ( (:fast_copy!, false), (:unsafe_fast_copy!, true) )
@eval begin

# x = y
# =====
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, incy::Int, n::Int)
    $isunsafe || begin
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    mul = max(abs(incx), abs(incy))
    if n < $NLIM_COPY1*mul || n*mul > $NLIM_COPY2
        # src: copy_foroop.jl
incx < 0 && (ix = ix+(n-1)*abs(incx))
incy < 0 && (iy = iy+(n-1)*abs(incy))
@inbounds for i = 0:n-1
    x[ix+i*incx] = y[iy+i*incy]
end
    else
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
    end
    return x
end
# inceq
function ($f)(x::Array{$elty}, ix::Int, incx::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inceq.jl
incy = incx
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    mul = abs(incx)
    if n < $NLIM_COPY1*mul || n*mul > $NLIM_COPY2
        # src: copy_foroop_inceq.jl
incx = abs(incx)
d = iy - ix
@inbounds for i = ix:incx:ix+(n-1)*incx
    x[i] = y[d+i]
end
    else
        # src: set2_inceq.jl
incy = incx
        # src: blas_copy.jl
px = convert(Ptr{$(elty)},x) + (ix-1)*sizeof($(elty))
py = convert(Ptr{$(elty)},y) + (iy-1)*sizeof($(elty))
ccall(($(string(fcopy)),libblas), Void,
    (Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}, Ptr{$(elty)}, Ptr{BlasInt}),
    &(n), py, &(incy), px, &(incx))
    end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, iy::Int, n::Int)
    $isunsafe || begin
        # src: set2_inc1.jl
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: set2_inc1.jl
incx = 1
incy = 1
    # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
    return x
end
# inc1ieq
function ($f)(x::Array{$elty}, ix::Int, y::Array{$elty}, n::Int)
    $isunsafe || begin
        # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
        # src: fast_check2.jl
(0 != incx && 0 != incy) || throw(ArgumentError("zero increment"))
(0 < ix && 0 < iy) || throw(BoundsError())
ix+(n-1)*abs(incx) <= length(x) || throw(BoundsError())
iy+(n-1)*abs(incy) <= length(y) || throw(BoundsError())
    end
    # src: set2_inc1ieq.jl
iy = ix
incx = 1
incy = 1
    # src: c_memcpy.jl
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
py = convert(Ptr{$(elty)},y) + (iy-1)*selty
ccall(:memcpy, Ptr{Void}, (Ptr{Void}, Ptr{Void}, Uint),
    px, py, n*selty)
    return x
end

end # eval begin
end # for

## FILL METHODS

for (f, isunsafe) in ( (:fast_fill!, false), (:unsafe_fast_fill!, true) )
@eval begin

# x = a
# =======
# general
function ($f)(x::Array{$elty}, ix::Int, incx::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    # src: fill_foroop.jl
incx = abs(incx)
@inbounds for i = ix:incx:ix-1+n*incx
    x[i] = a
end
    return x
end
# inc1
function ($f)(x::Array{$elty}, ix::Int, a::$elty, n::Int)
    $isunsafe || begin
        # src: set1_inc1.jl
incx = 1
        # src: fast_check1.jl
0 < incx || throw(ArgumentError("non-positive increment"))
0 < ix || throw(BoundsError())
ix+(n-1)*incx <= length(x) || throw(BoundsError())
    end
    if a == 0 && n > $NLIM_FILL && ZEROFLOAT
        # src: set1_inc1.jl
incx = 1
        # src: c_memset.jl
a::Int32 = a
selty = sizeof($(elty))
px = convert(Ptr{$(elty)},x) + (ix-1)*selty
ccall(:memset, Ptr{Void}, (Ptr{Void}, Int32, Csize_t),
    px, a, n*selty)
    else
        # src: fill_foroop_inc1.jl
@inbounds for i = ix:ix-1+n
    x[i] = a
end
    end
    return x
end

end # eval begin
end # for

end # for

end # module


text:
@testset "compiler_new" begin

eqs = [:(foo(0) = log(a(0))+b(0)/x(-1)), :(bar(0) = c(1)+u*d(1))]
args = [(:a, -1), (:a, 0), (:b, 0), (:c, 0), (:c, 1), (:d, 1)]
params = [:u]
defs = Dict(:x=>:(a(0)/(1-c(1))))
targets = [(:foo, 0), (:bar, 0)]
funname = :myfun

flat_args = [(:a, 0), (:b, 1), (:c, -1)]
grouped_args = OrderedDict(:x=>[(:a, -1),], :y=>[(:a, 0), (:b, 0), (:c, 0)], :z=>[(:c, 1), (:d, 1)])
flat_params = [:beta, :delta]
grouped_params = Dict(:p => [:u])

args2 = vcat(args, targets)::Vector{Tuple{Symbol,Int}}

ff = Dolang.FunctionFactory(eqs, args, params, targets=targets, defs=defs,
                            funname=funname)

ff2 = Dolang.FunctionFactory(eqs, grouped_args, params, targets=targets, defs=defs,
                            funname=funname)

# we're all static arrays
x = @SVector [1.0]
y = @SVector [2.0, 0.5, 0.3]
z = @SVector [0.4, 0.2]
p = @SVector [0.1]
#
# we're all vectors
xv = Vector(x)
yv = Vector(y)
zv = Vector(z)
pv = Vector(p)

# we're list of points
N = 5
x_vec = reshape(reinterpret(SVector{1,Float64}, 1 .+ rand(1,N)), (N,))
y_vec = reshape(reinterpret(SVector{3,Float64}, rand(3,N)), (N,))
z_vec = reshape(reinterpret(SVector{2,Float64}, rand(2,N)), (N,))
p_vec = reshape(reinterpret(SVector{1,Float64}, rand(1,N)), (N,))

# as matrices
x_mat = Dolang.from_SA(x_vec)
y_mat = Dolang.from_SA(y_vec)
z_mat = Dolang.from_SA(z_vec)
p_mat = Dolang.from_SA(p_vec)

@testset "Static Arrays Functions" begin

    @testset "Flat Function Factory constructors" begin
        @test begin
            Dolang.FlatFunctionFactory(ff)
            Dolang.FlatFunctionFactory(ff2)
            true
        end
    end

    @testset "Kernel" begin
        fff = Dolang.FlatFunctionFactory(ff2)
        cc = Dolang.gen_kernel(fff, [0, 1, 3], funname=:kernel)
        kernel = Core.eval(Dolang, cc)
        res = @inferred kernel(x, y, z, p)
        @test isa(res, Tuple{SVector{2,Float64},SMatrix{2,1,Float64},SMatrix{2,2,Float64}})
    end

    @testset "GuFun" begin
        fff = Dolang.FlatFunctionFactory(ff2)
        cc = Dolang.gen_gufun(fff, [0, 1, 3], funname=:gufun)
        gufun = Core.eval(Dolang, cc)
        # behaves like a kernel
        res = @inferred gufun(x, y, z, p)
        @test isa(res, Tuple{SVector{2,Float64},SMatrix{2,1,Float64},SMatrix{2,2,Float64}})
        # doesn't like allocated vectors but does it anyway
        res = @inferred gufun(xv, yv, zv, pv)
        @test isa(res, Tuple{Vector{Float64}, Matrix{Float64}, Matrix{Float64}})
        # vectorized calls
        res = @inferred gufun(x_vec, y_vec, z_vec, p_vec)
        @test isa(res,
            Tuple{Array{SVector{2,Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}}
        )
        @test (length(res[1]) == length(res[2]) == length(x_vec))
        res2 = @inferred gufun(x, y_vec, z, p_vec)
        @test isa(res2,
            Tuple{Array{SVector{2,Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}}
        )
        res3 = @inferred gufun(x_vec, y, z, p)
        @test isa(res3,
            Tuple{Array{SVector{2,Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}}
        )
        # vectorized calls with matrix inputs

        # TODO: the method below is not type-inferrable!
        # res4 = @inferred gufun(x_mat, y_mat, z_mat, p_mat)
        res4 = gufun(x_mat, y_mat, z_mat, p_mat)
        @test isa(res4, Tuple{Matrix{Float64},Array{Float64,3},Array{Float64,3}})

        # TODO: the method below is not type-inferrable!
        # res5 = @inferred gufun(x_mat, y_mat, z_mat, pv)
        res5 = gufun(x_mat, y_mat, z_mat, pv)
        @test isa(res5, Tuple{Matrix{Float64},Array{Float64,3},Array{Float64,3}})

        res6 = @inferred gufun(xv, yv, zv, pv)
        @test isa(res6, Tuple{Vector{Float64},Matrix{Float64},Matrix{Float64}})

    end

    @testset "Generated Gufuns" begin
        fff = Dolang.FlatFunctionFactory(ff2)
        cc = Dolang.gen_generated_gufun(fff; funname=:gengufun)
        gengufun = Core.eval(Dolang, cc)
        # behaves like a kernel
        res0 = @inferred gengufun(x, y, z, p)
        @test isa(res0, SVector{2,Float64})

        deriv_specs = [
            (Val{0}(), Val{1}(), Val{3}()),
            Val{(0, 1, 3)},
            Val{(0, 1, 3)}()
        ]
        for derivative in deriv_specs

            res = @inferred gengufun(derivative, x, y, z, p)
            @test isa(res, Tuple{SVector{2,Float64},SMatrix{2,1,Float64},SMatrix{2,2,Float64}})

            # doesn't like allocated vectors but does it anyway
            res = @inferred gengufun(derivative, xv, yv, zv, pv)
            @test isa(res, Tuple{Array{Float64,1},Array{Float64,2},Array{Float64,2}})

            # vectorized calls
            res = @inferred gengufun(derivative, x_vec, y_vec, z_vec, p_vec)
            @test isa(res, Tuple{Array{SVector{2,Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}})
            @test (length(res[1]) == length(res[2]) == length(x_vec))

            res2 = @inferred gengufun(derivative, x, y_vec, z, p_vec)
            @test isa(res2, Tuple{Array{SVector{2,Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}})

            res3 = @inferred gengufun(derivative, x_vec, y, z, p)
            @test isa(res3, Tuple{Array{SVector{2, Float64},1},Array{StaticArrays.SArray{Tuple{2,1},Float64,2,2},1},Array{StaticArrays.SArray{Tuple{2,2},Float64,2,4},1}})
        end

        @test isa(@inferred(gengufun(Val{0}, x, y, z, p)), SVector{2,Float64})
        @test isa(@inferred(gengufun(Val{1}, x, y, z, p)), StaticArrays.SArray{Tuple{2,1},Float64,2,2})
        @test isa(@inferred(gengufun(Val{(1,)}, x, y, z, p)), Tuple{StaticArrays.SArray{Tuple{2,1},Float64,2,2}})

    end

end
end

ext: jl

max_issues_repo_name: EconForge/Dolang


text:
struct AliasSampler{T<:AbstractFloat,SF<:Random.Sampler{T},SI<:Random.Sampler{Int}} <:
       Random.Sampler{Int}
    acceptance::Vector{T}
    aliastable::Vector{Int}
    floatsampler::SF
    integersampler::SI
end

function AliasSampler(
    prediction::AbstractVector{<:Real},
    floatsampler::Random.Sampler{T},
    integersampler::Random.Sampler{Int},
) where {T<:AbstractFloat}
    sum_prediction = sum(prediction)
    sum_prediction ≈ one(sum_prediction) || error("predictions have to be normalized")

    n = length(prediction)
    acceptance = Vector{T}(undef, n)
    aliastable = Vector{Int}(undef, n)
    StatsBase.make_alias_table!(prediction, sum_prediction, acceptance, aliastable)

    return AliasSampler(acceptance, aliastable, floatsampler, integersampler)
end

function Base.rand(rng::Random.AbstractRNG, sampler::AliasSampler)
    idx = rand(rng, sampler.integersampler)
    x = rand(rng, sampler.floatsampler)
    @inbounds r = x < sampler.acceptance[idx] ? idx : sampler.aliastable[idx]
    return r
end


text:
function 한자찾기(사전::한자사전, 한자::AbstractString)
  한자타입(사전[한자]...)
end

function 독음으로찾기(사전::한자사전, 찾을독음::AbstractString)
  결과 = []
  for (한자, (코드,정체,간체,독음,획수)) in 사전
    if occursin(찾을독음, 독음)
      push!(결과, 한자)
    end
  end
  결과
end

for (func, block) in [
    (:한글보기, (사전,글자) -> 한자찾기(사전, 글자).독음),
    (:정체보기, (사전,글자) -> begin
      한자 = 한자찾기(사전, 글자)
      isempty(한자.정체) ? 글자 : 한자.정체
    end),
    (:간체보기, (사전,글자) -> begin
      한자 = 한자찾기(사전, 글자)
      isempty(한자.간체) ? 글자 : 한자.간체
    end),
  ]
  @eval begin
    function $(func)(사전::한자사전, 문장::AbstractString)
      결과 = []
      for 글자 in split(문장, "")
        if haskey(사전, 글자)
          push!(결과, $block(사전, 글자))
        else
          push!(결과, 글자)
        end
      end
      join(결과)
    end
  end
end

for (func, block) in [
    (:한자옆에한글보기, (사전, 결과, 글자, 단어에추가) -> begin
      한자 = 한자찾기(사전, 글자)
      push!(결과, 글자)
      단어에추가(한자.독음)
    end),
    (:한글옆에한자보기, (사전, 결과, 글자, 단어에추가) -> begin
      한자 = 한자찾기(사전, 글자)
      push!(결과, 한자.독음)
      단어에추가(글자)
    end),
    (:한자옆에간체보기, (사전, 결과, 글자, 단어에추가) -> begin
      한자 = 한자찾기(사전, 글자)
      push!(결과, 글자)
      단어에추가(isempty(한자.간체) ? 글자 : 한자.간체)
    end),
    (:간체옆에정체보기, (사전, 결과, 글자, 단어에추가) -> begin
      한자 = 한자찾기(사전, 글자)
      push!(결과, isempty(한자.간체) ? 글자 : 한자.간체)
      단어에추가(isempty(한자.정체) ? 글자 : 한자.정체)
    end),
  ]
  @eval begin
    function $(func)(사전::한자사전, 문장::AbstractString)
      결과 = []
      찾음 = false
      단어 = []
      function 단어에추가(대상)
        push!(단어, 대상)
        찾음 = true
      end
      function 찾은단어넣기()
        if 찾음
          push!(결과, "($(join(단어)))")
          찾음 = false
          단어 = []
        end
      end
      for 글자 in split(문장, "")
        if haskey(사전, 글자)
          $block(사전, 결과, 글자, 단어에추가)
        else
          찾은단어넣기()
          push!(결과, 글자)
        end
      end
      찾은단어넣기()
      join(결과)
    end
  end
end


text:
import Base.Threads.@spawn
@everywhere using IterTools
N = 5
op = [:+, :-, :*, :/, :^]
fun = [:sqrt, :factorial]

@everywhere function create_expr(elem)
    if isa(elem, Tuple)
        return Expr(:call, [create_expr(e) for e in elem]...)
    else
        return elem
    end
end
#create_expr(t::Tuple) = Expr(:call, [if typeof(x) == Tuple create_expr(x) else x end for x in t]...)

@everywhere function perform(all_values)
    push!(all_values, N)
    push!(all_values, 10*N + N)
    ll = length(all_values)
    size_op = length(op) * ll^2
    size_fun = (ll + size_op) * length(fun)
    #size_fun = size_fun + size_op*length(op)* length(fun)
    size_fun = ll * length(fun)

    all_op = Iterators.product(op, all_values, all_values)
    fun_val = Iterators.product(fun, all_values)
    #fun_op = Iterators.product(fun, all_op)
    #op_fun = Iterators.product(op, all_op, fun)
    #all_op_fun = Iterators.flatten((all_op, fun_val, fun_op, op_fun))
    op_fun = Iterators.flatten((all_op, fun_val))

    all_exprs = Vector{Expr}(undef, size_op + size_fun)

    for (i, t) in enumerate(op_fun)
        all_exprs[i] = create_expr(t)
    end

    return all_exprs
end

@everywhere function expands(expr::Expr)
    args = []
    for a in expr.args
        if typeof(a) == Expr
            args = [args..., expands(a)...]
        else
            push!(args, a)
        end
    end
    return args
end

@everywhere count_val(val, expr::Expr) = count(i -> i == val, expands(expr))
@everywhere count_N(expr::Expr) = count_val(N, expr) + 2*count_val(N*10 + N, expr)

@everywhere lessThanN(expr::Expr) = count_N(expr) <= N

@everywhere equal_N(expr::Expr) = count_N(expr) == N

@everywhere function factorial_div_sqrt(expr::Expr)
    if expr.args[1] != :factorial
        return true
    else
        expr_vals = expands(expr)
        return !(any(f -> f == :/, expr_vals) || any(f -> f == :sqrt, expr_vals))
    end
end

@everywhere function pow_negative_big(expr::Expr)
    if expr.args[1] != :^
        return true
    else
        return (eval(expr.args[2]) >= 0) && (0 <= eval(expr.args[3]) <= 16)
    end
end

@everywhere function sqrt_fact_negative(expr::Expr)
    if !in(expr.args[1], [:sqrt, :factorial])
        return true
    else
        expr_vals = expands(expr)
        if any(f -> in(f, [:/, :^, :sqrt]), expr_vals)
            return false
        end
        return eval(expr.args[2]) >= 0
    end
end

@everywhere function big_factorials(expr::Expr)
    if expr.args[1] != :factorial
        return true
    else
        expr_vals = expands(expr)
        if any(f -> in(f, [:/, :^, :sqrt]), expr_vals)
            return false
        end
        return eval(expr.args[2]) < 20
    end
end

@everywhere function reduce_equivalents(expr_vals)
    dd = Dict()
    for (value, expr) in expr_vals
        if haskey(dd, value)
            push!(dd[value], expr)
        else
            dd[value] = [expr]
        end
    end

    new_expr_vals = []
    for (val, li_expr) in dd
        nb_n_done = []
        for expr in li_expr
            nb_n = count_N(expr)
            if !in(nb_n, nb_n_done)
                push!(new_expr_vals, (val, expr))
                push!(nb_n_done, nb_n)
            end
        end
    end

    return new_expr_vals
end

@everywhere no_pow_mul(expr::Expr) = !in(expr.args[2], [:*, :^])
@everywhere no_sqrt(expr::Expr) = !in(expr.args[1], [:sqrt])
@everywhere two_add(expr::Expr) = count(i -> i == :+, expands(expr)) == 2
@everywhere no_55(expr::Expr) = count(i -> i == 55, expands(expr)) == 0

function main(number)
    all_values = [] #[N, 10*N + N]
    for i = 1:N
        all_values = perform(all_values)

        all_values = Iterators.filter(lessThanN, all_values)
        all_values = Iterators.filter(factorial_div_sqrt, all_values)
        all_values = Iterators.filter(big_factorials, all_values)
        all_values = Iterators.filter(sqrt_fact_negative, all_values)
        all_values = Iterators.filter(pow_negative_big, all_values)

        expr_vals = map(expr -> (eval(expr), expr), all_values)
        expr_vals = filter(item -> - N*100 < item[1] < N*100, expr_vals)
        expr_vals = filter(item -> floor(item[1]) == item[1], expr_vals)
        expr_vals = reduce_equivalents(expr_vals)

        all_exprs = []
        for (value, expr) in expr_vals
            if eval(expr) == number && equal_N(expr)
                return expr
            end
            push!(all_exprs, expr)
        end
        all_values = all_exprs # map(item -> item[2], expr_vals)

    end
end

for i = 60:100
    expr = main(i)
    println(eval(expr),": ", expr)
end

max_stars_repo_name: atthom/CheckFourFours

max_issues_repo_name: atthom/CheckFourFours

max_forks_repo_name: atthom/CheckFourFours


text:
#Shorthand for JavaCall.iterator(iterable_collection).
macro iter(collection)
   return :( JavaCall.iterator($collection) )
end

#Convert to string and then println.
macro print_string(object)
   return :( println(to_string($object)) )
end

#Convert to list and then println.
macro print_list(object)
   return :( println(to_list($object)) )
end

#Convert to Penn and then println.
macro print_penn(object)
   return :( println(pennstring($object)) )
end


text:
using Test, YaoBlocks

@testset "test constant gates" begin
    include("const_gate.jl")
end

@testset "test identity gates" begin
    include("identity_gate.jl")
end

@testset "test phase gate" begin
    include("phase_gate.jl")
end

@testset "test shift gate" begin
    include("shift_gate.jl")
end

@testset "test rotation gate" begin
    include("rotation_gate.jl")
end

@testset "test time evolution" begin
    include("time_evolution.jl")
end

@testset "test general matrix gate" begin
    include("general_matrix_gate.jl")
end

@testset "test measure" begin
    include("measure.jl")
end

# it does nothing
@test chsubblocks(X, Y) === X


text:
function get_cache_slots_number(::Type{OrekitConfiguration})
    return jcall(OrekitConfiguration, "getCacheSlotsNumber", jint, ())
end

function set_cache_slots_number(::Type{OrekitConfiguration}, arg0::jint)
    return jcall(OrekitConfiguration, "setCacheSlotsNumber", void, (jint,), arg0)
end


text:
"""
Crude but more-or-less effective thing that prints out a struct (which may contain other structs) as
a markdown table.
"""
function to_md_table( f; exclude=[], depth=0 ) :: String
    F = typeof(f)
    @assert isstructtype( F )
    names = fieldnames(F)
    prinames = []
    structnames = []

    for n in names
        v = getfield(f,n)
        if n in exclude
            ;
        elseif isstructtype( typeof(v))
            push!(structnames, n )
        else
            push!(prinames, n )
        end
    end
    s = """
    |            |              |
    |:-----------|-------------:|
    """
    for n in prinames
        v = getfield(f,n)
        s *= "|$n|$v|\n"
    end
    s *= "\n\n"
    depth += 1
    for n in structnames
        v = getfield(f,n)
        s *= "#"*repeat( "#", depth ) * "$n\n"
        s *= to_md_table( v, exclude=exclude, depth=depth )
    end
    return s;
end

max_stars_repo_path: src/function to_string( thing ).jl

max_stars_repo_name: grahamstark/ScottishTaxBenefitModel

max_issues_repo_path: src/function to_string( thing ).jl

max_issues_repo_name: grahamstark/ScottishTaxBenefitModel

max_forks_repo_path: src/function to_string( thing ).jl

max_forks_repo_name: grahamstark/ScottishTaxBenefitModel


text: # placeholder for prepare_test_data.jl


text:
function fma(a::T, b::T, c::T) where {F<:IEEEFloat, T<:DoubleFloat{F}}
    hi, lo = fma_2(a, b, c)
    return DoubleFloat(hi, lo)
end

function fma_4(a::T, b::T, c::T) where {F<:IEEEFloat, T<:DoubleFloat{F}}
    c_hi = HI(c) * onehalf(F)
    c_lo = LO(c) * onehalf(F)

    fma1 = fma(HI(a), HI(b), c_hi)
    fma2 = fma(HI(a), LO(b), c_hi)
    fma3 = fma(LO(a), HI(b), c_lo)
    fma4 = fma(LO(a), LO(b), c_lo)

    return add_4(fma1, fma2, fma3, fma4)
end

function fma_(a::T, b::T, c::T) where {F<:IEEEFloat, T<:DoubleFloat{F}}
    c_hi = HI(c) * half(F)
    c_lo = LO(c) * half(F)

    fma1 = fma(HI(a), HI(b), c_hi)
    fma2 = fma(HI(a), LO(b), c_hi)
    fma3 = fma(LO(a), HI(b), c_lo)
    fma4 = fma(LO(a), LO(b), c_lo)

    return add_3(fma1, fma2, fma3, fma4)
end

fma_3(a::T, b::T, c::T) where {F<:IEEEFloat, T<:DoubleFloat{F}} =
    fma_(a, b, c)

function fma_2(a::T, b::T, c::T) where {F<:IEEEFloat, T<:DoubleFloat{F}}
    c_hi = halfx(HI(c))
    c_lo = halfx(LO(c))

    fma1 = fma(HI(a), HI(b), c_hi)
    fma2 = fma(HI(a), LO(b), c_hi)
    fma3 = fma(LO(a), HI(b), c_lo)
    fma4 = fma(LO(a), LO(b), c_lo)

    return add_2(fma1, fma2, fma3, fma4)
end


text:
using Documenter, RemoteREPL

makedocs(;
    modules=[RemoteREPL],
    format=Documenter.HTML(),
    pages=[
        "Overview"   => "index.md",
        "Tutorial"   => "tutorial.md",
        "How To"     => "howto.md",
        "Reference"  => "reference.md",
        "Design"     => "design.md",
    ],
    repo="https://github.com/c42f/RemoteREPL.jl/blob/{commit}{path}#L{line}",
    sitename="RemoteREPL.jl",
    authors = "Chris Foster and contributors: https://github.com/c42f/RemoteREPL.jl/graphs/contributors"
)

deploydocs(;
    repo="github.com/c42f/RemoteREPL.jl",
    devbranch="main",
    push_preview=true
)


text:
# Julia wrapper for header: cpl_error.h
# Automatically generated using Clang.jl

"""
    CPLEmergencyError(const char * pszMessage) -> void

Fatal error when things are bad.

### Parameters
* **pszMessage**: the error message to report.
"""
function cplemergencyerror(arg1)
    ccall((:CPLEmergencyError, libgdal), Cvoid, (Cstring,), arg1)
end

"""
    CPLErrorReset() -> void

Erase any traces of previous errors.
"""
function cplerrorreset()
    ccall((:CPLErrorReset, libgdal), Cvoid, ())
end

"""
    CPLGetLastErrorNo() -> CPLErrorNum

Fetch the last error number.

### Returns
the error number of the last error to occur, or CPLE_None (0) if there are no posted errors.
"""
function cplgetlasterrorno()
    ccall((:CPLGetLastErrorNo, libgdal), CPLErrorNum, ())
end

"""
    CPLGetLastErrorType() -> CPLErr

Fetch the last error type.

### Returns
the error type of the last error to occur, or CE_None (0) if there are no posted errors.
"""
function cplgetlasterrortype()
    ccall((:CPLGetLastErrorType, libgdal), CPLErr, ())
end

"""
    CPLGetLastErrorMsg() -> const char *

Get the last error message.

### Returns
the last error message, or NULL if there is no posted error message.
"""
function cplgetlasterrormsg()
    unsafe_string(ccall((:CPLGetLastErrorMsg, libgdal), Cstring, ()))
end

"""
    CPLGetErrorCounter() -> GUInt32

Get the error counter.

### Returns
the error counter.
"""
function cplgeterrorcounter()
    ccall((:CPLGetErrorCounter, libgdal), GUInt32, ())
end

"""
    CPLGetErrorHandlerUserData(void) -> void *

Fetch the user data for the error context.

### Returns
the user data pointer for the error context
"""
function cplgeterrorhandleruserdata()
    ccall((:CPLGetErrorHandlerUserData, libgdal), Ptr{Cvoid}, ())
end

"""
    CPLErrorSetState(CPLErr eErrClass,
                     CPLErrorNum err_no,
                     const char * pszMsg) -> void

Restore an error state, without emitting an error.
"""
function cplerrorsetstate(eErrClass, err_no, pszMsg)
    ccall((:CPLErrorSetState, libgdal), Cvoid, (CPLErr, CPLErrorNum, Cstring), eErrClass, err_no, pszMsg)
end

"""
    CPLCleanupErrorMutex() -> void
"""
function cplcleanuperrormutex()
    ccall((:CPLCleanupErrorMutex, libgdal), Cvoid, ())
end

"""
    CPLLoggingErrorHandler(CPLErr eErrClass,
                           CPLErrorNum nError,
                           const char * pszErrorMsg) -> void

Error handler that logs into the file defined by the CPL_LOG configuration option, or stderr otherwise.
"""
function cplloggingerrorhandler(arg1, arg2, arg3)
    ccall((:CPLLoggingErrorHandler, libgdal), Cvoid, (CPLErr, CPLErrorNum, Cstring), arg1, arg2, arg3)
end

"""
    CPLDefaultErrorHandler(CPLErr eErrClass,
                           CPLErrorNum nError,
                           const char * pszErrorMsg) -> void

Default error handler.
"""
function cpldefaulterrorhandler(arg1, arg2, arg3)
    ccall((:CPLDefaultErrorHandler, libgdal), Cvoid, (CPLErr, CPLErrorNum, Cstring), arg1, arg2, arg3)
end

"""
    CPLQuietErrorHandler(CPLErr eErrClass,
                         CPLErrorNum nError,
                         const char * pszErrorMsg) -> void

Error handler that does not do anything, except for debug messages.
"""
function cplquieterrorhandler(arg1, arg2, arg3)
    ccall((:CPLQuietErrorHandler, libgdal), Cvoid, (CPLErr, CPLErrorNum, Cstring), arg1, arg2, arg3)
end

"""
    CPLTurnFailureIntoWarning(int bOn) -> void

Whether failures should be turned into warnings.
"""
function cplturnfailureintowarning(bOn)
    aftercare(ccall((:CPLTurnFailureIntoWarning, libgdal), Cvoid, (Cint,), bOn))
end

"""
    CPLSetErrorHandler(CPLErrorHandler pfnErrorHandlerNew) -> CPLErrorHandler

Install custom error handler.

### Parameters
* **pfnErrorHandlerNew**: new error handler function.

### Returns
returns the previously installed error handler.
"""
function cplseterrorhandler(arg1)
    ccall((:CPLSetErrorHandler, libgdal), CPLErrorHandler, (CPLErrorHandler,), arg1)
end

"""
    CPLSetErrorHandlerEx(CPLErrorHandler pfnErrorHandlerNew,
                         void * pUserData) -> CPLErrorHandler

Install custom error handle with user's data.

### Parameters
* **pfnErrorHandlerNew**: new error handler function.
* **pUserData**: User data to carry along with the error context.

### Returns
returns the previously installed error handler.
"""
function cplseterrorhandlerex(arg1, arg2)
    ccall((:CPLSetErrorHandlerEx, libgdal), CPLErrorHandler, (CPLErrorHandler, Ptr{Cvoid}), arg1, arg2)
end

"""
    CPLPushErrorHandler(CPLErrorHandler pfnErrorHandlerNew) -> void

Push a new CPLError handler.

### Parameters
* **pfnErrorHandlerNew**: new error handler function.
"""
function cplpusherrorhandler(arg1)
    ccall((:CPLPushErrorHandler, libgdal), Cvoid, (CPLErrorHandler,), arg1)
end

"""
    CPLPushErrorHandlerEx(CPLErrorHandler pfnErrorHandlerNew,
                          void * pUserData) -> void

Push a new CPLError handler with user data on the error context.

### Parameters
* **pfnErrorHandlerNew**: new error handler function.
* **pUserData**: User data to put on the error context.
"""
function cplpusherrorhandlerex(arg1, arg2)
    ccall((:CPLPushErrorHandlerEx, libgdal), Cvoid, (CPLErrorHandler, Ptr{Cvoid}), arg1, arg2)
end

"""
    CPLSetCurrentErrorHandlerCatchDebug(int bCatchDebug) -> void

Set if the current error handler should intercept debug messages, or if they should be processed by the previous handler.

### Parameters
* **bCatchDebug**: FALSE if the current error handler should not intercept debug messages
"""
function cplsetcurrenterrorhandlercatchdebug(bCatchDebug)
    ccall((:CPLSetCurrentErrorHandlerCatchDebug, libgdal), Cvoid, (Cint,), bCatchDebug)
end

"""
    CPLPopErrorHandler() -> void

Pop error handler off stack.
"""
function cplpoperrorhandler()
    ccall((:CPLPopErrorHandler, libgdal), Cvoid, ())
end

"""
    _CPLAssert(const char * pszExpression,
               const char * pszFile,
               int iLine) -> void

Report failure of a logical assertion.
"""
function _cplassert(arg1, arg2, arg3)
    aftercare(ccall((:_CPLAssert, libgdal), Cvoid, (Cstring, Cstring, Cint), arg1, arg2, arg3))
end


text:
import Base: ==

"""
    IsEqual(fn)

Takes a funciton of the form `fn(x)::Bool`
"""
struct IsEqual
    fn::Function
end

(==)(x::IsEqual, y) = x.fn(y)
(==)(y , x::IsEqual) = x == y

issymbol(x) = x isa Symbol

"""
    get_symbol(expr)

Return `x` for expressions of form `x::Type` otherwise return nothing
"""
function get_symbol(expr)
    if expr == Expr(:(::), IsEqual(issymbol), IsEqual(x->true))
        expr.args[1]
    else
        nothing
    end
end

"""
    get_type(x)

Return `T` if an expresion is of the form `:(x::Type{T})` or `:(::Type{T})` when `T` is a symbol otherwise returns `nothing`
"""
function get_type(expr)
    if expr == Expr(:(::), Expr(:curly, :Type, IsEqual(issymbol)))
        expr.args[1].args[2]
    elseif expr == Expr(:(::), IsEqual(issymbol) , Expr(:curly, :Type, IsEqual(issymbol)))
        expr.args[2].args[2]
    else
        nothing
    end
end

"""
    getargs_dottilde(x)

Return the arguments `L` and `R`, if `x` is an expression of the form `L .~ R` or
`(~).(L, R)`, or `nothing` otherwise.
"""
function getargs_dottilde(expr)
    any_arg = IsEqual(x->true)
    # Check if the expression is of the form `L .~ R`.
    if expr == Expr(:call, :.~, any_arg, any_arg)
         expr.args[2], expr.args[3]
    # Check if the expression is of the form `(~).(L, R)`.
    elseif expr == Expr(:., :~, Expr(:tuple, any_arg, any_arg))
        expr.args[2].args[1], expr.args[2].args[2]
    else
        nothing
    end
end

"""
    getargs_tilde(x)

Return the arguments `L` and `R`, if `x` is an expression of the form `L ~ R`, or `nothing`
otherwise.
"""
function getargs_tilde(expr)
    any_arg = IsEqual(x->true)
    if expr == Expr(:call, :~, any_arg, any_arg)
         expr.args[2], expr.args[3]
    else
        nothing
    end
end

############################################
# Julia 1.2 temporary fix - Julia PR 33303 #
############################################
function to_namedtuple_expr(syms, vals=syms)
    if length(syms) == 0
        nt = :(NamedTuple())
    else
        nt_type = Expr(:curly, :NamedTuple,
            Expr(:tuple, QuoteNode.(syms)...),
            Expr(:curly, :Tuple, [:(Core.Typeof($x)) for x in vals]...)
        )
        nt = Expr(:call, :($(DynamicPPL.namedtuple)), nt_type, Expr(:tuple, vals...))
    end
    return nt
end

if VERSION == v"1.2"
    @eval function namedtuple(::Type{NamedTuple{names, T}}, args::Tuple) where {names, T <: Tuple}
        if length(args) != length(names)
            throw(ArgumentError("Wrong number of arguments to named tuple constructor."))
        end
        # Note T(args) might not return something of type T; e.g.
        # Tuple{Type{Float64}}((Float64,)) returns a Tuple{DataType}
        $(Expr(:splatnew, :(NamedTuple{names,T}), :(T(args))))
    end
else
    function namedtuple(::Type{NamedTuple{names, T}}, args::Tuple) where {names, T <: Tuple}
        return NamedTuple{names, T}(args)
    end
end

#####################################################
# Helper functions for vectorize/reconstruct values #
#####################################################

vectorize(d::UnivariateDistribution, r::Real) = [r]
vectorize(d::MultivariateDistribution, r::AbstractVector{<:Real}) = copy(r)
vectorize(d::MatrixDistribution, r::AbstractMatrix{<:Real}) = copy(vec(r))

# NOTE:
# We cannot use reconstruct{T} because val is always Vector{Real} then T will be Real.
# However here we would like the result to be specifric type, e.g. Array{Dual{4,Float64}, 2},
# otherwise we will have error for MatrixDistribution.
# Note this is not the case for MultivariateDistribution so I guess this might be lack of
# support for some types related to matrices (like PDMat).
reconstruct(d::UnivariateDistribution, val::AbstractVector) = val[1]
reconstruct(d::MultivariateDistribution, val::AbstractVector) = copy(val)
function reconstruct(d::MatrixDistribution, val::AbstractVector)
    return reshape(copy(val), size(d))
end
function reconstruct!(r, d::Distribution, val::AbstractVector)
    return reconstruct!(r, d, val)
end
function reconstruct!(r, d::MultivariateDistribution, val::AbstractVector)
    r .= val
    return r
end
function reconstruct(d::Distribution, val::AbstractVector, n::Int)
    return reconstruct(d, val, n)
end
function reconstruct(d::UnivariateDistribution, val::AbstractVector, n::Int)
    return copy(val)
end
function reconstruct(d::MultivariateDistribution, val::AbstractVector, n::Int)
    return copy(reshape(val, size(d)[1], n))
end
function reconstruct(d::MatrixDistribution, val::AbstractVector, n::Int)
    tmp = reshape(val, size(d)[1], size(d)[2], n)
    orig = [tmp[:, :, i] for i in 1:size(tmp, 3)]
    return orig
end
function reconstruct!(r, d::Distribution, val::AbstractVector, n::Int)
    return reconstruct!(r, d, val, n)
end
function reconstruct!(r, d::MultivariateDistribution, val::AbstractVector, n::Int)
    r .= val
    return r
end

# Uniform random numbers with range 4 for robust initializations
# Reference: https://mc-stan.org/docs/2_19/reference-manual/initialization.html
randrealuni() = 4 * rand() - 2
randrealuni(args...) = 4 .* rand(args...) .- 2

const Transformable = Union{PositiveDistribution,UnitDistribution,TransformDistribution,
                            SimplexDistribution,PDMatDistribution}
istransformable(dist) = false
istransformable(::Transformable) = true

#################################
# Single-sample initialisations #
#################################

inittrans(dist::UnivariateDistribution) = invlink(dist, randrealuni())
inittrans(dist::MultivariateDistribution) = invlink(dist, randrealuni(size(dist)[1]))
inittrans(dist::MatrixDistribution) = invlink(dist, randrealuni(size(dist)...))

################################
# Multi-sample initialisations #
################################

inittrans(dist::UnivariateDistribution, n::Int) = invlink(dist, randrealuni(n))
function inittrans(dist::MultivariateDistribution, n::Int)
    return invlink(dist, randrealuni(size(dist)[1], n))
end
function inittrans(dist::MatrixDistribution, n::Int)
    return invlink(dist, [randrealuni(size(dist)...) for _ in 1:n])
end


text:
"""
    greens_lorentzian(κ)

Evaluate the Green's function corresponding to a Lorentzian lineshape at `κ = (t-τ)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> greens_lorentzian((t-τ)/T2s)
4.5399929762484854e-5
```
"""
function greens_lorentzian(κ)
    exp(-κ)
end

"""
    dG_o_dT2s_x_T2s_lorentzian(κ)

Evaluate the derivative of Green's function, corresponding to a Lorentzian lineshape, wrt. `T2s` at `κ = (t-τ)/T2s` and multiply it by `T2s`.

The multiplication is added so that the function merely depends on `κ = (t-τ)/T2s`. The actual derivative is given by `dG_o_dT2s_x_T2s_lorentzian((t-τ)/T2s)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> dGdT2s = dG_o_dT2s_x_T2s_lorentzian((t-τ)/T2s)/T2s
45.39992976248485
```
"""
function dG_o_dT2s_x_T2s_lorentzian(κ)
    greens_lorentzian(κ) * κ
end

"""
    greens_gaussian(κ)

Evaluate the Green's function corresponding to a Gaussian lineshape at `κ = (t-τ)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> greens_gaussian((t-τ)/T2s)
1.9287498479639178e-22
```
"""
function greens_gaussian(κ)
    exp(-κ^2 / 2)
end

"""
    dG_o_dT2s_x_T2s_gaussian(κ)

Evaluate the derivative of Green's function, corresponding to a Gaussian lineshape, wrt. `T2s` at `κ = (t-τ)/T2s` and multiply it by `T2s`.

The multiplication is added so that the function merely depends on `κ = (t-τ)/T2s`. The actual derivative is given by `dG_o_dT2s_x_T2s_gaussian((t-τ)/T2s)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> dGdT2s = dG_o_dT2s_x_T2s_gaussian((t-τ)/T2s)/T2s
1.9287498479639177e-15
```
"""
function dG_o_dT2s_x_T2s_gaussian(κ)
    greens_gaussian(κ) * κ^2
end

"""
    greens_superlorentzian(κ)

Evaluate the Green's function corresponding to a super-Lorentzian lineshape at `κ = (t-τ)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> greens_superlorentzian((t-τ)/T2s)
0.1471246868094442
```
"""
function greens_superlorentzian(κ)
    quadgk(ζ -> exp(- κ^2 * (3ζ^2 - 1)^2 / 8), 0, sqrt(1/3), 1, order = 100)[1]
end

"""
    dG_o_dT2s_x_T2s_superlorentzian(κ)

Evaluate the derivative of Green's function, corresponding to a super-Lorentzian lineshape, wrt. `T2s` at `κ = (t-τ)/T2s` and multiply it by `T2s`.

The multiplication is added so that the function merely depends on `κ = (t-τ)/T2s`. The actual derivative is given by `dG_o_dT2s_x_T2s_superlorentzian((t-τ)/T2s)/T2s`.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> dGdT2s = dG_o_dT2s_x_T2s_superlorentzian((t-τ)/T2s)/T2s
15253.095033670965
```
"""
function dG_o_dT2s_x_T2s_superlorentzian(κ)
    quadgk(ζ -> exp(-κ^2 * (3ζ^2 - 1)^2 / 8) * (κ^2 * (3ζ^2 - 1)^2 / 4), 0, sqrt(1/3), 1, order = 100)[1]
end

"""
    interpolate_greens_function(f, κmin, κmax)

Interpolate the Green's function f in the range between κmin and κmax.

The interpolation uses the ApproxFun.jl package that incorporates Chebyshev polynomials and ensures an approximation to machine precision.

# Examples
```jldoctest
julia> t = 100e-6;

julia> τ = 0;

julia> T2s = 10e-6;

julia> greens_superlorentzian((t-τ)/T2s)
0.1471246868094442

julia> Gint = interpolate_greens_function(greens_superlorentzian, 0, 20);

julia> Gint((t-τ)/T2s)
0.14712468680944407
```
"""
function interpolate_greens_function(f, κmin, κmax)
    Fun(f, κmin..κmax)
    # x = κmin : 0.001 : κmax
    # return CubicSplineInterpolation(x, f.(x))
end


text:
mutable struct Points
    lonlats::Vector{Tuple{Float64,Float64}}
    placeids::Vector{String}
    originalindex::Vector{Int}
end

function getpoints(points::Vector, originalindex::Vector{Int})
    for (i,p) in enumerate(points)
        if haskey(p, "originalIndex")
            originalindex[p["originalIndex"]+1] = i
        end
    end
    Points(
        [(l["location"]["longitude"],l["location"]["latitude"]) for l in points],
        [l["placeId"] for l in points],
        originalindex
    )
end

function matchroads(snap::Bool, npoints::Int; kwargs...)
    resultjson = if snap
        requestjson("snapToRoads", "https://roads.googleapis.com/v1/", ""; kwargs...)
    else
        requestjson("nearestRoads", "https://roads.googleapis.com/v1/", ""; kwargs...)
    end
    getpoints(resultjson["snappedPoints"], fill(0, npoints))
end

"""

The Google Maps Roads API takes takes up to 100 independent coordinates, and
returns the closest road segment for each point. The points passed do not need
to be part of a continuous path.

If you are working with sequential GPS points, set `snap=true`.

## Parameters
* `key`: your google maps API key
* `lonlats`: the list of points to match

## Keyword Arguments
* `snap`: if you are working with sequential GPS points
* `interpolate`: (only for `snap=true`) Whether to include all points forming
    the full road geometry. When true, additional interpolated points will also
    be returned, resulting in a path that smoothly follows the geometry of the
    road, even around corners and through tunnels. Interpolated paths will most
    likely contain more points than the original path. Defaults to `false`.

## Remarks
The snapping algorithm works best for points that are not too far apart. If you
observe odd snapping behavior, try creating paths that have points closer
together. To ensure the best snap-to-road quality, you should aim to provide
paths on which consecutive pairs of points are within 300m of each other. This
will also help in handling any isolated, long jumps between consecutive points
caused by GPS signal loss, or noise.
"""
function matchroads(
        key::String,
        lonlats::Vector;
        snap::Bool = false,
        maxcoords::Int = 100,
        kwargs...
    )
    npoints = length(lonlats)
    if npoints > 100
        warn("The API takes up to 100 lon,lats. Ignoring coords 101...$npoints")
        npoints = min(maxcoords, npoints)
    end
    points = join([string(lonlats[i][2], ",", lonlats[i][1]) for i in 1:npoints], "|")
    if snap
        matchroads(snap, npoints, key = key, path = points; kwargs...)
    else
        matchroads(snap, npoints, key = key, points = points; kwargs...)
    end
end


text:
#
# Keep the books for aasol
#
mutable struct ItStatsA{T<:Real}
    condhist::Array{T,1}
    alphanorm::Array{T,1}
    history::Array{T,1}
end

function ItStatsA(rnorm)
    ItStatsA([1.0],[1.0],[rnorm])
end

#
# Collect the stats at the end of the iteration
#
function CollectStats(ItData::ItStatsA)
    stats = (condhist=ItData.condhist[2:end],
                    alphanorm=ItData.alphanorm[2:end])
return stats
end

function updateStats!(ItData::ItStatsA, condhist, alphanorm)
      append!(ItData.condhist, condhist)
      append!(ItData.alphanorm, alphanorm)
end

function updateHist!(ItData::ItStatsA, rnorm)
      append!(ItData.history, rnorm)
end

#
# Initialize Anderson iteration
#
function Anderson_Init(x0, Vstore, m, maxit, beta, keepsolhist)
blocksize=1024
(0.0 < abs(beta) <= 1) || error("abs(beta) must be in (0,1]")
sol=copy(x0)
n=length(x0)
(mv, nv) = size(Vstore)
mv == n || error("Vstore needs ", n, " rows")
(nv >= 2 * (m + 1)) || error("Vstore needs ", 2 * m + 4, " columns")
#
# Just in case you are reusing Vstore for several problems, I will
# reinitialize it to zero.
#
Vstore .= 0.0
if m==0
   Qd=[]; QP=[]; DG=[];
   nvblock=1
else
QP = @views Vstore[:, 1:m]
DG = @views Vstore[:, m+1:2*m]
if (nv >= 3*m+3)
    (Qd = @views Vstore[:,2*m+1:3*m-1])
    nvblock=3*m
else
    @warn "Low storage mode"
    Qd=zeros(blocksize,m-1)
    nvblock=2*m+1
end
end
gx = Anderson_vector_Init(Vstore,nvblock)
df = Anderson_vector_Init(Vstore,nvblock+1)
dg = Anderson_vector_Init(Vstore,nvblock+2)
res = Anderson_vector_Init(Vstore,nvblock+3)
keepsolhist ? (solhist = solhistinit(n, maxit, sol)) : (solhist = [])
return(sol, gx, df, dg, res, DG, QP, Qd, solhist)
end

function Anderson_vector_Init(Vstore,nvblock)
gx = @views Vstore[:,nvblock]
return gx
end

#
# Figure out what idid and errcode are. Boring but must be done.
# This is a lot simpler than for Newton. There are no linesearches
# or Krylov iterations to keep track of.
#
function AndersonOK(resnorm, tol, k, toosoon)
idid=(resnorm <= tol)
idid ? (errcode = 0) : (errcode = 10)
toosoon && (errcode = -1)
idid || println("Failure to converge in aasol.jl")
toosoon && println("Iteration terminates on entry to aasol.jl")
return (idid, errcode)
end

"""
falpha(alpha,theta,mk)
Map thetas to alphas for stats
"""
function falpha(alpha, theta, mk)
    alpha[1] = theta[1]
    for ia = 2:mk
        alpha[ia] = theta[ia] - theta[ia-1]
    end
    alpha[mk+1] = 1.0 - theta[mk]
    return norm(alpha, 1)
end


text:
using DiffusionDefinition
using Test
using StaticArrays, Random
DD = DiffusionDefinition

@testset "DiffusionDefinition.jl" begin
    @diffusion_process TestDiffusion begin
        :dimensions
        process --> 3
        wiener --> 5

        :parameters
        _ --> (3, Float64)
        param --> Float64
        stem --> (2, Int64)
        (theta, alpha) --> (2, Int32)
        (beta, gamma) --> Float32
        (yota, zeta) --> (Float64, Int64)

        :aux_info
        T --> Float64
        vT --> SArray{Tuple{4},Float64,1,4}
        #=
        :conjugate
        phi(t, u) --> (
            (p1*u[2]*(u[2]-1),u[1],u[1]),
            (p3*u[5]*t,u[2],u[4]),
        )
        nonhypo(x) --> x
        num_non_hypo --> 3
        =#
        :additional
        statespace --> LowerBoundedStateSpace((1,2,3), (0.0, 1.0, 0.0))
        constdiff --> false
        linear --> false
        eltype --> Float64
    end a b c d

    param_tuple = (1.0, 2.0, 3.0, 4.0, 5, 6, Int32(7), Int32(8),
                      Float32(9.0), Float32(10.0), 11.0, 12)
    end_point_info_tuple = (1.0, (@SVector [1.0, 2.0, 3.0, 4.0]))
    P = TestDiffusion(param_tuple..., end_point_info_tuple...)

    @test DD.nonhypo(:x, P) == :x
    @testset "satisfying bounds" begin
        @test !DD.bound_satisfied(P, [0.1, 2.0, 0.0])
        @test DD.bound_satisfied(P, [0.1, 2.0, 0.1])
        @test !DD.bound_satisfied(P, [0.1, 1.0, 0.1])
    end
    @test !DD.constdiff(P)
    @test typeof(P) <: DD.DiffusionProcess
    @test !(typeof(P) <: DD.LinearDiffusion)
    @test eltype(P) == Float64

    param_names = (:p1, :p2, :p3, :param, :stem1, :stem2, :theta, :alpha, :beta, :gamma, :yota, :zeta, :T, :vT)
    @test DD.parameter_names(P) == param_names
    #TODO introduce this back ?
    #@test DD.parameter_names(typeof(P)) == (:p1, :p2, :p3, :param, :stem1, :stem2, :theta, :alpha, :beta, :gamma, :yota, :zeta)
    params_is = DD.parameters(P)
    params_should = Dict([n=>p for (n,p) in zip(param_names, param_tuple)]...)
    for p in keys(params_should)
        @test params_is[p] == params_should[p]
    end
    @test DD.end_point_info_names(P) == (:T, :vT)
    epi_is = DD.end_point_info(P)
    epi_should = Dict(
        :T=>end_point_info_tuple[1],
        :vT=>end_point_info_tuple[2],
    )
    for p in keys(epi_should)
        @test epi_is[p] == epi_should[p]
    end

    @diffusion_process TestDiffusion2 begin
        :dimensions
        process --> 3
        wiener --> 5

        :parameters
        _ --> (3, Float64)

        :aux_info
        T --> Float64
        vT --> SArray{Tuple{4},Float64,1,4}

    end

    param_tuple = (1.0, 2.0, 3.0)
    end_point_info_tuple = (1.0, (@SVector [1.0, 2.0, 3.0, 4.0]))

    P = TestDiffusion2(param_tuple..., end_point_info_tuple...)

    @test DD.parameter_names(P) == (:p1, :p2, :p3, :T, :vT)

    params_is = DD.parameters(P)
    params_should = Dict(:p1 => 1.0, :p2 => 2.0, :p3 => 3.0, :T => 1.0, :vT => (@SVector [1.0, 2.0, 3.0, 4.0]))
    for p in keys(params_should)
        @test params_is[p] == params_should[p]
    end
    @test DD.end_point_info_names(P) == (:T, :vT)
    epi_is = DD.end_point_info(P)
    epi_should = Dict(
        :T=>end_point_info_tuple[1],
        :vT=>end_point_info_tuple[2],
    )
    for p in keys(epi_should)
        @test epi_is[p] == epi_should[p]
    end
    P1 = TestDiffusion2(1.0, 1.0, 1.0, end_point_info_tuple...)
    P2 = TestDiffusion2(1.0, 1.0, 1.0, 1.0, (@SVector [1.0, 1.0, 1.0, 1.0]))

    #@test DD.update_params(P, [1.0, 1.0, 1.0]...) == P1

    #@test DD.clone(P, [1.0, 1.0, 1.0]..., [1.0, (@SVector [1.0, 1.0, 1.0, 1.0]), (@SVector [1.0, 1.0, 1.0])]...) == P2
end

@testset "State space restrictions" begin
    lb = LowerBoundedStateSpace((1,3), (-2.0, 3.0))
    ub = UpperBoundedStateSpace((1,3), (4.0, 10.0))
    b = BoundedStateSpace(lb, ub)
    nb = UnboundedStateSpace()

    a = [-1,0,1.0,4.0]
    v = [-1,0,1.0,20.0]
    @test DD._bound_satisfied(nb, a)
    @test !DD._bound_satisfied(lb, a)
    @test DD._bound_satisfied(ub, a)
    @test !DD._bound_satisfied(b, a)
    @test DD._bound_satisfied(nb, v)
    @test !DD._bound_satisfied(lb, v)
    @test DD._bound_satisfied(ub, v)
    @test !DD._bound_satisfied(b, v)
end

@testset "utility functions" begin
    @test DD.remove_curly(Vector{Float64}) == Array
    @test DD.get_curly(Vector{Float64}) == (:Float64, 1)
end

@testset "Buffers" begin
    d, m = 4,5
    el = Float64
    sb = DD.StandardEulerBuffer(zeros(el,d), zeros(el,(d,m)), zeros(el,m))
    #@test size(sb) == (d + d + d*m + m,)
    @test length(sb.b) == d
    @test size(sb.σ) == (d, m)
    @test length(sb.dW) == m
    #@test eltype(sb) == el
    #sb2 = similar(sb)
    #@test sb2 !== sb
    #sb3 = similar(sb, SVector{3,Float64})
    #@test size(sb3) == (d + d*m + m,)
    #@test sb3.data == zeros(SVector{3,Float64}, d+d*m+m)
    #=
    @testset "indexing StandardEulerBuffer" begin
        sb.b[1] = 1.0
        #@test sb[:] == sb.data
        sb[2:5] .= [2.0, 3.0, 4.0, 5.0]
        sb[[6,8]] .= [6.0, 8.0]
        sb[15] = 15.0
        sb[25] = 25.0

        @test sb.b == [1.0, 2.0, 3.0, 4.0]
        @test sb.y == [5.0, 6.0, 0.0, 8.0]
        @test sb.σ == [
            0.0 0.0 0.0 0.0 25.0;
            0.0 0.0 0.0 0.0 0.0;
            0.0 15.0 0.0 0.0 0.0;
            0.0 0.0 0.0 0.0 0.0
        ]
        @test sb.dW == [0.0, 0.0, 0.0, 0.0, 0.0]
    end
    =#
    lb = DD.LinearDiffBuffer(zeros(el,d), zeros(el,(d,m)), zeros(el,m), zeros(el,(d,d)))
    #@test size(lb) == (d + d + d*m + d*d + m,)
    @test length(lb.b) == d
    @test size(lb.B) == (d,d)
    @test size(lb.σ) == (d, m)
    @test length(lb.dW) == m
end

@testset "setting parameters and cloning" begin
    @load_diffusion LotkaVolterra
    @load_diffusion LotkaVolterraAux

    θ = [2.0/3.0, 4.0/3.0, 1.0, 1.0, 0.2, 0.3]

    P_target = LotkaVolterra(θ...)
    P_aux = LotkaVolterraAux(θ..., 0.0, 1.0, (@SVector [1.0, 2.0]))

    p_names = (:α, :β, :γ, :δ, :σ1, :σ2)
    pars = Dict([pn=>p for (pn,p) in zip(p_names,θ)])
    ep_names = (:t0, :T, :vT)

    @test DD.parameter_names(LotkaVolterra) == p_names
    @test DD.parameter_names(P_target) == p_names
    @test DD.parameters(P_target) == pars

    @test DD.const_parameter_names(LotkaVolterra) == tuple()
    @test DD.const_parameter_names(P_target) == tuple()
    @test DD.const_parameters(P_target) == Dict{Any,Any}()

    @test DD.var_parameter_names(LotkaVolterra) == p_names
    @test DD.var_parameter_names(P_target) == p_names
    @test DD.var_parameters(P_target) == pars

    @test DD.end_point_info_names(LotkaVolterra) == tuple()
    @test DD.end_point_info_names(P_target) == tuple()
    @test DD.end_point_info(P_target) == Dict{Any,Any}()

    @test DD.parameter_names(LotkaVolterraAux) == Tuple(vcat(p_names..., ep_names...))
    @test DD.parameter_names(P_aux) == Tuple(vcat(p_names..., ep_names...))

    #TODO finish below...
    DD.parameters(P_aux)

    DD.const_parameter_names(LotkaVolterraAux)
    DD.const_parameter_names(P_aux)
    DD.const_parameters(P_aux)

    DD.var_parameter_names(LotkaVolterraAux)
    DD.var_parameter_names(P_aux)
    DD.var_parameters(P_aux)

    DD.end_point_info_names(LotkaVolterraAux)
    DD.end_point_info_names(P_aux)
    DD.end_point_info(P_aux)

    DD.set_parameters!(P_target, [10.0, 20.0, 30.0, 40.0], [2=>:δ, 4=>:α])
    @test P_target.δ == 20.0
    @test P_target.α == 40.0
end

@testset "Sampling trajectories" begin

    Random.seed!(3)
    N, dt, T = 30, 0.0001, 10.0
    tt = collect(0.0:dt:T)
    quad_var(x) = sum(abs2.(diff(x)))

    tr = trajectory(tt, SVector{N,Float64})
    wr = Wiener()
    rand!(wr, tr)

    tr_mut = trajectory(tt, Vector{Float64}, N)
    wr_mut = Wiener(N)
    rand!(wr_mut, tr_mut)

    @test maximum([abs(quad_var(map(x->x[i], tr.x))/T-1.0) for i in 1:N])<0.05
    @test maximum([abs(quad_var(map(x->x[i], tr_mut.x))/T-1.0) for i in 1:N])<0.05

    # an example of a pre-defined diffusion
    @load_diffusion LotkaVolterraAux

    α, β, γ, δ, σ1, σ2 = 2.0/3.0, 4.0/3.0, 1.0, 1.0, 0.2, 0.3
    lv_aux = LotkaVolterraAux(
        α, β, γ, δ, σ1, σ2,
        0.0, 1.0, zero(DD.ℝ{2})
    )
    N, dt, T = 2, 0.0001, 1000.0
    tt = collect(0.0:dt:T)

    WW = trajectory(tt, SVector{N,Float64})
    wr = Wiener()
    rand!(wr, WW)

    v0 = @SVector [1.0, 2.0]
    XX = trajectory(tt, SVector{N,Float64})
    DD.solve!(XX, WW, lv_aux, v0)

    mean(x) = sum(x)/length(x)
    x1 = map(x->x[1], XX.x)
    x2 = map(x->x[2], XX.x)

    @test abs(mean(x1)-γ/δ) < 0.05
    @test abs(mean(x2)-α/β) < 0.05
    @test abs(quad_var(x1)/T-σ1^2) < 0.01
    @test abs(quad_var(x2)/T-σ2^2) < 0.01

    # and in-place simulation
    function DD.B!(buffer, t, P::LotkaVolterraAux)
        #B = DD.get_B(buffer)
        buffer.B[1,1] = 0.0
        buffer.B[1,2] = -P.β*P.γ/P.δ
        buffer.B[2,1] = P.α*P.δ/P.β
        buffer.B[2,2] = 0.0
    end

    function DD.β!(buffer, t, P::LotkaVolterraAux)
        #b = DD.get_b(buffer)
        buffer.b[1] = P.γ/P.δ*P.α
        buffer.b[2] = -P.α/P.β*P.γ
    end

    function DD.σ!(buffer, t, P::LotkaVolterraAux)
        #σ = DD.get_σ(buffer)
        buffer.σ[1,1] = P.σ1
        buffer.σ[1,2] = 0.0
        buffer.σ[2,1] = 0.0
        buffer.σ[2,2] = P.σ2
    end

    WW = trajectory(tt, Vector{Float64}, N)
    wr = Wiener(N)
    rand!(wr, WW)

    XX = trajectory(tt, Vector{Float64}, N)
    v0 = [1.0, 2.0]
    buffer = DD.LinearDiffBuffer{Vector{Float64}}(lv_aux)

    DD.solve!(XX, WW, lv_aux, v0, buffer)

    x1 = map(x->x[1], XX.x)
    x2 = map(x->x[2], XX.x)
    @test abs(mean(x1)-γ/δ) < 0.05
    @test abs(mean(x2)-α/β) < 0.05
    @test abs(quad_var(x1)/T-σ1^2) < 0.05
    @test abs(quad_var(x2)/T-σ2^2) < 0.05
end


text:
## 2.1 - Strings

# Q1
hi1 = repeat("hi", 1000)
hi2 = "hi"^1000

@assert hi1 == hi2

# Q2

a = 3
b = 4

str1 = "$a + $b"
str2 = "$(a+b)"

## 2.2 - Data Structures

# Q1

a_ray = [1, 2, 3]
push!(a_ray, 4)
println(a_ray)
pop!(a_ray)
println(a_ray)

# Q2

# This dictionary can only contain strings (i.e. type is Dict{String, String})
myphonebook = Dict("Jenny" => "867-5309", "Ghostbusters" => "555-2368")

# This dictionary can contain any type of object (i.e. type is Dict{String, Any})
flexible_phonebook = Dict("Jenny" => 8675309, "Ghostbusters" => "555-2368")

max_stars_repo_path: Exercise Solutions/2 - Basics.jl

max_issues_repo_path: Exercise Solutions/2 - Basics.jl

max_forks_repo_path: Exercise Solutions/2 - Basics.jl


text:
# This file was generated, do not modify it. # hide
pipe = Pipeline(
    coercer = X -> coerce(X, :age=>Continuous),
    one_hot_encoder = OneHotEncoder(),
    transformed_target_model = TransformedTargetModel(
        model = KNNRegressor(K=3);
        target=UnivariateStandardizer()
    )
)

max_issues_repo_name: alan-turing-institute/MLJTutorials

max_forks_repo_name: alan-turing-institute/MLJTutorials


text:
module Klimakoffer

using LinearAlgebra: \, lu
using SparseArrays: sparse
using UnPack: @unpack
using Requires: @require

export Mesh, Model, Discretization, compute_equilibrium!, compute_evolution!, set_co2_concentration!

function __init__()
    # Enable features that depend on the availability of the Plots package
    @require Plots="91a5bcdd-55d7-5caf-9e0b-520d859cae80" begin
      using .Plots: Plots
      include("visualization.jl")
    end

    @require Makie="ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" begin
      using .Makie: Makie
    end
end

# Include additional files
include("mesh.jl")
include("model.jl")
include("discretization.jl")
include("numerics.jl")

end # module


text:
#-------------------------------------------------------------------
#* EMSO Model Library (EML) Copyright (C) 2004 - 2007 ALSOC.
#*
#* This LIBRARY is free software; you can distribute it and/or modify
#* it under the therms of the ALSOC FREE LICENSE as available at
#* http://www.enq.ufrgs.br/alsoc.
#*
#* EMSO Copyright (C) 2004 - 2007 ALSOC, original code
#* from http://www.rps.eng.br Copyright (C) 2002-2004.
#* All rights reserved.
#*
#* EMSO is distributed under the therms of the ALSOC LICENSE as
#* available at http://www.enq.ufrgs.br/alsoc.
#*----------------------------------------------------------------------
#* Author: Gerson Balbueno Bicca
#* $Id$
#*--------------------------------------------------------------------
type heater
	heater()=begin
		PP=outers.PP
		NComp=outers.NComp
		new(
			DanaPlugin(Dict{Symbol,Any}(
				:Brief=>"Physical Properties",
				:Type=>"PP"
			)),
			DanaInteger(Dict{Symbol,Any}(
				:Brief=>"Number of Components"
			)),
			DanaSwitcher(Dict{Symbol,Any}(
				:Brief=>"Option for Display Phase Equilibrium K-values",
				:Valid=>["yes","no"],
				:Default=>"yes"
			)),
			power(Dict{Symbol,Any}(
				:Brief=>"Actual Duty",
				:Symbol=>"Q_{Duty}"
			)),
			fraction(Dict{Symbol,Any}(
				:Brief=>"Vapor fraction Outlet Stream",
				:Symbol=>"V_{frac}"
			)),
			fraction(Dict{Symbol,Any}(
				:Brief=>"Liquid fraction Outlet Stream",
				:Symbol=>"L_{frac}"
			)),
			fill(DanaReal(Dict{Symbol,Any}(
				:Brief=>"Phase Equilibrium K-values",
				:Lower=>1E-30,
				:Upper=>1E30,
				:Symbol=>"K_{value}"
			)),(NComp)),
			positive(Dict{Symbol,Any}(
				:Brief=>"Pressure Ratio",
				:Symbol=>"P_{ratio}"
			)),
			press_delta(Dict{Symbol,Any}(
				:Brief=>"Pressure Drop",
				:DisplayUnit=>"kPa",
				:Symbol=>"\\Delta P"
			)),
			stream(Dict{Symbol,Any}(
				:Brief=>"Inlet Stream",
				:PosX=>0,
				:PosY=>0.40,
				:Symbol=>"^{in}"
			)),
			streamPH(Dict{Symbol,Any}(
				:Brief=>"Outlet Stream",
				:PosX=>1,
				:PosY=>0.40,
				:Symbol=>"^{out}"
			)),
			power(Dict{Symbol,Any}(
				:Brief=>"Inlet Heat Stream",
				:PosX=>0.5,
				:PosY=>1,
				:Symbol=>"_{out}"
			)),
			[
				:(Outlet.F = Inlet.F),
				:(Outlet.F*Outlet.z([1 : NComp]) = Inlet.F*Inlet.z([1 : NComp])),
				:(Vfrac = Outlet.v),
				:(Lfrac = 1-Vfrac),
				:(Duty = Outlet.F*Outlet.h - Inlet.F*Inlet.h),
				:(Outlet.P = Inlet.P - Pdrop),
				:(Outlet.P = Inlet.P * Pratio),
				:(Kvalue*Outlet.x = Outlet.y),
				:(Kvalue = 1),
				:(Heat = Duty),
			],
			[
				"Flow","Composition","Vapor fraction Outlet Stream","Liquid fraction Outlet Stream","Heat Duty","Pressure Drop","Pressure Ratio","K-values Phase Equilibrium","K-values Phase Equilibrium","Duty Specification",
			],
			[:PP,:NComp,:Kvalues,],
			[:Duty,:Vfrac,:Lfrac,:Kvalue,:Pratio,:Pdrop,:Inlet,:Outlet,:Heat,]
		)
	end
	PP::DanaPlugin
	NComp::DanaInteger
	Kvalues::DanaSwitcher
	Duty::power
	Vfrac::fraction
	Lfrac::fraction
	Kvalue::Array{DanaReal}
	Pratio::positive
	Pdrop::press_delta
	Inlet::stream
	Outlet::streamPH
	Heat::power
	equations::Array{Expr,1}
	equationNames::Array{String,1}
	parameters::Array{Symbol,1}
	variables::Array{Symbol,1}
	attributes::Dict{Symbol,Any}
end
export heater
function setEquationFlow(in::heater)
	addEquation(1)
	addEquation(2)
	addEquation(3)
	addEquation(4)
	addEquation(5)
	addEquation(6)
	addEquation(7)
	let switch=Kvalues
		# Fix for better convergence !!!
		if switch=="yes"
			addEquation(8)
		elseif switch=="no"
			addEquation(9)
		end
	end
	addEquation(10)
end
function atributes(in::heater,_::Dict{Symbol,Any})
	fields::Dict{Symbol,Any}=Dict{Symbol,Any}()
	fields[:Pallete]=true
	fields[:Icon]="icon/heater"
	fields[:Brief]="Heater"
	fields[:Info]="Determines thermal and phase conditions of an outlet stream.

== Specify ==
* The Inlet stream
* Specify:
**The outlet temperature and the outlet pressure or
** The outlet temperature and the inlet energy stream or
** The outlet pressure and the inlet energy stream
"
	drive!(fields,_)
	return fields
end
heater(_::Dict{Symbol,Any})=begin
	newModel=heater()
	newModel.attributes=atributes(newModel,_)
	newModel
end


text:
function find_triple(target)
    a_limit = Int(floor(target/(2+sqrt(2))))
    for a in 1:a_limit
        b = ((target^2/2)-target*a)/(target-a)
        if isinteger(b)
            println(Int(prod([a, b, sqrt(a^2+b^2)])))
        end
    end
end

find_triple(1000)

max_stars_repo_name: HarrisonGreen/Project-Euler-Solutions

max_issues_repo_name: HarrisonGreen/Project-Euler-Solutions

max_forks_repo_name: HarrisonGreen/Project-Euler-Solutions


text:
prompt(msg, yes::Bool = false) = prompt(stdin, msg, yes)

function prompt(io::IO, msg, yes::Bool = false)
    print(msg)

    if yes
        println(" Yes.")
    else
        print(" [Y/n] ")
        run(`stty raw`)
        input = read(io, Char)
        run(`stty cooked`)
        println()
        input in ['Y', 'y', '\n', '\r'] || return false
    end

    return true
end

sysimg() = "libcomonicon.$(Libdl.dlext)"
sysimg(name) = "lib$name.$(Libdl.dlext)"

"""
    default_exename()
Default Julia executable name: `joinpath(Sys.BINDIR, Base.julia_exename())`
"""
default_exename() = joinpath(Sys.BINDIR, Base.julia_exename())

dot_julia() = first(DEPOT_PATH)

"""
    default_julia_bin()
Return the default path to `.julia/bin`.
"""
default_julia_bin() = joinpath(first(DEPOT_PATH), "bin")

"""
    default_julia_fpath()
Return the default path to `.julia/completions`
"""
default_julia_fpath() = joinpath(first(DEPOT_PATH), "completions")