Skip to content

Instantly share code, notes, and snippets.

@ericphanson
Created August 3, 2019 12:47
Show Gist options
  • Star 1 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save ericphanson/21774ecdecc1788715858466de2ff4a2 to your computer and use it in GitHub Desktop.
Save ericphanson/21774ecdecc1788715858466de2ff4a2 to your computer and use it in GitHub Desktop.
Various attempts to differentiate through `eigmin` in Julia

Attempts to differentiate eigmin

Julia 1.1, with Zygote master

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteMaster/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2 #master (https://github.com/FluxML/Zygote.jl.git)
  [700de1a5] ZygoteRules v0.1.0 #master (https://github.com/FluxML/ZygoteRules.jl)

julia> versioninfo()
Julia Version 1.1.1
Commit 55e36cc (2019-05-16 04:10 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin15.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -1.1949325143627607
eigmin_vec(v) = -1.9065242723547482
Zygote.gradient(eigmin_vec, v) = ([0.247035, 0.00872636, 0.208667, -0.284921, 0.00872636, 0.000678498, 0.00532077, -0.0194134, 0.208667, 0.00532077, 0.187612, -0.188899, -0.284921, -0.0194134, -0.188899, 0.564674],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)```

Results:
```julia
ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2339#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2339#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((ArnoldiMethod.#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((getfield(ArnoldiMethod, Symbol("#kw##partialschur"))())))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:555
 [25] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [26] include(::Module, ::String) at ./sysimg.jl:29
 [27] exec_options(::Base.JLOptions) at ./client.jl:267
 [28] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -0.536820849533713
test_eigmin(A) = -0.5368208495337133

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Applications/Julia-1.1.app/Contents/Resources/julia/share/julia/stdlib/v1.1/LinearAlgebra/src/blas.jl:456
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:87
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [5] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53
 [6] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [7] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [8] _schur! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [9] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [10] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [11] _eigvals! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [12] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:6 [inlined]
 [13] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [14] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:31
 [15] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:37
 [16] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:46
 [17] top-level scope at show.jl:555
 [18] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [19] include(::Module, ::String) at ./sysimg.jl:29
 [20] exec_options(::Base.JLOptions) at ./client.jl:267
 [21] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -1.1634638294015776
test_eigmin(A) = -1.1634638294015784

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:598 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:608
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2339#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2339#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((getfield(IterativeSolvers, Symbol("#kw##lobpcg!"))())))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((IterativeSolvers.#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((getfield(IterativeSolvers, Symbol("#kw##lobpcg"))())))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((IterativeSolvers.#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((IterativeSolvers.lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((IterativeSolvers.#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((IterativeSolvers.lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:555
 [34] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [35] include(::Module, ::String) at ./sysimg.jl:29
 [36] exec_options(::Base.JLOptions) at ./client.jl:267
 [37] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -0.7722779425593689
test_eigmin(A) = -0.7722779425591891

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2339#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2339#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((KrylovKit.#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((KrylovKit.#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:555
 [30] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [31] include(::Module, ::String) at ./sysimg.jl:29
 [32] exec_options(::Base.JLOptions) at ./client.jl:267
 [33] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -0.6884119101624921
test_eigmin(A) = -0.6884119101624923

Julia 1.2, with Zygote master

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteMaster/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2 #master (https://github.com/FluxML/Zygote.jl.git)
  [700de1a5] ZygoteRules v0.1.0 #master (https://github.com/FluxML/ZygoteRules.jl)

julia> versioninfo()
Julia Version 1.2.0-rc2.0
Commit 9248bf7687 (2019-07-08 19:42 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin18.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -0.5164694205909903
eigmin_vec(v) = -1.6777098633424616
Zygote.gradient(eigmin_vec, v) = ([0.1213750771214421, 0.09754706520978354, -0.1535129476181751, -0.1298258062521984, 0.09754706520978354, 0.1470479791665443, 0.019031225931315898, -0.18848902588803446, -0.1535129476181751, 0.019031225931315898, 0.4895633862291685, -0.010356521541228042, -0.1298258062521984, -0.18848902588803446, -0.010356521541228042, 0.24201355748284503],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((#partialschur)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:576
 [25] include at ./boot.jl:328 [inlined]
 [26] include_relative(::Module, ::String) at ./loading.jl:1094
 [27] include(::Module, ::String) at ./Base.jl:31
 [28] exec_options(::Base.JLOptions) at ./client.jl:295
 [29] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -0.9015796088863606
test_eigmin(A) = -0.9015796088863606

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/blas.jl:456
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:87
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [5] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53
 [6] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [7] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [8] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [9] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [10] test_eigmin at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [11] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [12] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:31
 [13] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:37
 [14] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:46
 [15] top-level scope at show.jl:576
 [16] include at ./boot.jl:328 [inlined]
 [17] include_relative(::Module, ::String) at ./loading.jl:1094
 [18] include(::Module, ::String) at ./Base.jl:31
 [19] exec_options(::Base.JLOptions) at ./client.jl:295
 [20] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -1.0239568533992993
test_eigmin(A) = -1.0239568533992989

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:598 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:608
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((#lobpcg!)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((#lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:576
 [34] include at ./boot.jl:328 [inlined]
 [35] include_relative(::Module, ::String) at ./loading.jl:1094
 [36] include(::Module, ::String) at ./Base.jl:31
 [37] exec_options(::Base.JLOptions) at ./client.jl:295
 [38] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -1.0651622515340569
test_eigmin(A) = -1.0651622515338184

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:576
 [30] include at ./boot.jl:328 [inlined]
 [31] include_relative(::Module, ::String) at ./loading.jl:1094
 [32] include(::Module, ::String) at ./Base.jl:31
 [33] exec_options(::Base.JLOptions) at ./client.jl:295
 [34] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -1.1550796185683476
test_eigmin(A) = -1.1550796185683472

Julia 1.3, with Zygote master

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteMaster/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2 #master (https://github.com/FluxML/Zygote.jl.git)
  [700de1a5] ZygoteRules v0.1.0 #master (https://github.com/FluxML/ZygoteRules.jl)

julia> versioninfo()
Julia Version 1.3.0-alpha.0
Commit 6c11e7c2c4 (2019-07-23 01:46 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin18.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -1.4952733109588594
eigmin_vec(v) = -2.4145434500709424
Zygote.gradient(eigmin_vec, v) = ([0.1583789527794141, -0.0768770614777082, 0.20183444499552697, -0.21924769054434584, -0.0768770614777082, 0.0660300112281631, -0.11516350825885954, 0.029780006880861624, 0.20183444499552697, -0.11516350825885954, 0.26750794918760523, -0.23351242718403944, -0.21924769054434584, 0.029780006880861624, -0.23351242718403944, 0.508083086804818],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:612
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:625
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((#partialschur)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:582
 [25] include at ./boot.jl:328 [inlined]
 [26] include_relative(::Module, ::String) at ./loading.jl:1094
 [27] include(::Module, ::String) at ./Base.jl:31
 [28] exec_options(::Base.JLOptions) at ./client.jl:295
 [29] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -1.2817055585069217
test_eigmin(A) = -1.2817055585069215

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.3/LinearAlgebra/src/blas.jl:455
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:91
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [5] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53
 [6] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [7] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [8] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [9] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [10] test_eigmin at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [11] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [12] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:31
 [13] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:37
 [14] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:46
 [15] top-level scope at show.jl:582
 [16] include at ./boot.jl:328 [inlined]
 [17] include_relative(::Module, ::String) at ./loading.jl:1094
 [18] include(::Module, ::String) at ./Base.jl:31
 [19] exec_options(::Base.JLOptions) at ./client.jl:295
 [20] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -0.5321502950284704
test_eigmin(A) = -0.5321502950284704

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:612 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:622
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((#lobpcg!)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((#lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:582
 [34] include at ./boot.jl:328 [inlined]
 [35] include_relative(::Module, ::String) at ./loading.jl:1094
 [36] include(::Module, ::String) at ./Base.jl:31
 [37] exec_options(::Base.JLOptions) at ./client.jl:295
 [38] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -1.330893212569147
test_eigmin(A) = -1.3308932125684043

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:612
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:625
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/emit.jl:117
 [13] #s2317#1293 at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:20 [inlined]
 [14] #s2317#1293(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/SnR3X/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:582
 [30] include at ./boot.jl:328 [inlined]
 [31] include_relative(::Module, ::String) at ./loading.jl:1094
 [32] include(::Module, ::String) at ./Base.jl:31
 [33] exec_options(::Base.JLOptions) at ./client.jl:295
 [34] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -0.7056615828260215
test_eigmin(A) = -0.7056615828260214

Julia 1.1, with Zygote release

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteRelease/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2

julia> versioninfo()
Julia Version 1.1.1
Commit 55e36cc (2019-05-16 04:10 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin15.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -0.8599881111735685
eigmin_vec(v) = -1.8822909708572109
Zygote.gradient(eigmin_vec, v) = ([0.151512, 0.119427, -0.0856678, -0.246522, 0.119427, 0.343148, 0.0494105, -0.209114, -0.0856678, 0.0494105, 0.103352, 0.132439, -0.246522, -0.209114, 0.132439, 0.401987],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s79#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s79#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((ArnoldiMethod.#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((getfield(ArnoldiMethod, Symbol("#kw##partialschur"))())))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:555
 [25] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [26] include(::Module, ::String) at ./sysimg.jl:29
 [27] exec_options(::Base.JLOptions) at ./client.jl:267
 [28] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -1.355963664123503
test_eigmin(A) = -1.3559636641235033

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Applications/Julia-1.1.app/Contents/Resources/julia/share/julia/stdlib/v1.1/LinearAlgebra/src/blas.jl:456
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:87
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [5] _hessenberg! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53 [inlined]
 [6] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [7] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [8] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [9] _schur! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [10] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [11] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [12] _eigvals! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [13] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:6 [inlined]
 [14] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [15] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:31
 [16] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:37
 [17] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:46
 [18] top-level scope at show.jl:555
 [19] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [20] include(::Module, ::String) at ./sysimg.jl:29
 [21] exec_options(::Base.JLOptions) at ./client.jl:267
 [22] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -0.7091240336554688
test_eigmin(A) = -0.709124033655469

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:598 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:608
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s79#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s79#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((getfield(IterativeSolvers, Symbol("#kw##lobpcg!"))())))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((IterativeSolvers.#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((getfield(IterativeSolvers, Symbol("#kw##lobpcg"))())))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((IterativeSolvers.#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((IterativeSolvers.lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((IterativeSolvers.#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((IterativeSolvers.lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:555
 [34] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [35] include(::Module, ::String) at ./sysimg.jl:29
 [36] exec_options(::Base.JLOptions) at ./client.jl:267
 [37] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -0.7148727950599743
test_eigmin(A) = -0.7148727950594805

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:166 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:813 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s79#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s79#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((KrylovKit.#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((KrylovKit.#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((KrylovKit.eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:555
 [30] include_relative(::Module, ::String) at /Applications/Julia-1.1.app/Contents/Resources/julia/lib/julia/sys.dylib:?
 [31] include(::Module, ::String) at ./sysimg.jl:29
 [32] exec_options(::Base.JLOptions) at ./client.jl:267
 [33] _start() at ./client.jl:436
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -1.2316564609801013
test_eigmin(A) = -1.2316564609801015

Julia 1.2, with Zygote release

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteRelease/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2

julia> versioninfo()
Julia Version 1.2.0-rc2.0
Commit 9248bf7687 (2019-07-08 19:42 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin18.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -0.6430164197959672
eigmin_vec(v) = -1.6173855466810616
Zygote.gradient(eigmin_vec, v) = ([0.007739762987139389, 0.0421347696898412, -0.007472221526669799, -0.03931639904567614, 0.0421347696898412, 0.23227275719121965, -0.07881618049481061, -0.2019559426943517, -0.007472221526669799, -0.07881618049481061, 0.5098409572313541, -0.12124736047015983, -0.03931639904567614, -0.2019559426943517, -0.12124736047015983, 0.2501465225902866],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((#partialschur)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:576
 [25] include at ./boot.jl:328 [inlined]
 [26] include_relative(::Module, ::String) at ./loading.jl:1094
 [27] include(::Module, ::String) at ./Base.jl:31
 [28] exec_options(::Base.JLOptions) at ./client.jl:295
 [29] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -0.7709562710269048
test_eigmin(A) = -0.7709562710269056

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.2/LinearAlgebra/src/blas.jl:456
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:87
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [5] _hessenberg! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53 [inlined]
 [6] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [7] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [8] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [9] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [10] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [11] test_eigmin at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [12] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [13] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:31
 [14] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:37
 [15] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:46
 [16] top-level scope at show.jl:576
 [17] include at ./boot.jl:328 [inlined]
 [18] include_relative(::Module, ::String) at ./loading.jl:1094
 [19] include(::Module, ::String) at ./Base.jl:31
 [20] exec_options(::Base.JLOptions) at ./client.jl:295
 [21] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -0.7934792323543227
test_eigmin(A) = -0.7934792323543218

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:598 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:608
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((#lobpcg!)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((#lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:576
 [34] include at ./boot.jl:328 [inlined]
 [35] include_relative(::Module, ::String) at ./loading.jl:1094
 [36] include(::Module, ::String) at ./Base.jl:31
 [37] exec_options(::Base.JLOptions) at ./client.jl:295
 [38] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -1.3123829115843535
test_eigmin(A) = -1.3123829112135936

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:598
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:611
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] Type at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:576
 [30] include at ./boot.jl:328 [inlined]
 [31] include_relative(::Module, ::String) at ./loading.jl:1094
 [32] include(::Module, ::String) at ./Base.jl:31
 [33] exec_options(::Base.JLOptions) at ./client.jl:295
 [34] _start() at ./client.jl:464
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -1.2325376752252737
test_eigmin(A) = -1.2325376752252755

Julia 1.3, with Zygote release

Version information:

julia> Pkg.status()
    Status `~/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/ZygoteRelease/Project.toml`
  [ec485272] ArnoldiMethod v0.0.4
  [14197337] GenericLinearAlgebra v0.1.0
  [42fd0dbc] IterativeSolvers v0.8.1
  [0b1a1467] KrylovKit v0.3.4
  [e88e6eb3] Zygote v0.3.2

julia> versioninfo()
Julia Version 1.3.0-alpha.0
Commit 6c11e7c2c4 (2019-07-23 01:46 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin18.6.0)
  CPU: Intel(R) Core(TM) i5-5257U CPU @ 2.70GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, broadwell)

Method @adjoint (real rep of Hermitian).jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

# helper for `herm`
function f(x, y, i, j)
    if i == j
        return x + 0im
    elseif i > j
        return x + im * y
    elseif i < j
        return x - im * y
    end
end

# Create a complex `d` by `d` Hermitian matrix from a `d^2`-dimensional real vector `v`
function herm(v::AbstractVector)
    d = isqrt(length(v))
    @assert d^2 == length(v)
    m = reshape(v, d, d)
    m = [ f(m[i,j], m[j,i], i, j) for i = 1:d, j = 1:d ]
    return m
end

# Recover the vector representation
invherm(A::AbstractMatrix) = vec(real(A))

# Find the minimal eigenvalue using KrylovKit (from a vector representation)
function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigsolve(A, 1, :SR)[1][1] |> real
end

# Give an adjoint for `Zygote`, which technically only works in the non-degenerate case
Zygote.@adjoint function eigmin_vec(v::AbstractVector)
    A = herm(v)
    eigs = eigsolve(A, 1, :SR)
    eval = eigs[1][1]
    evec = eigs[2][1]
    adj = (evec * evec') |> invherm
    return eval |> real, c->(c * adj,)
end


A = rand(4,4); A = A + A';
v = invherm(A)

@show eigmin(A)
@show eigmin_vec(v)
@show Zygote.gradient(eigmin_vec, v)

Results:

eigmin(A) = -1.5506329143065933
eigmin_vec(v) = -2.7985341943349127
Zygote.gradient(eigmin_vec, v) = ([0.0505178879344079, 0.07913540401608901, -0.01936272471334829, -0.15323703912883263, 0.07913540401608901, 0.292160742770081, 0.13019912358366506, -0.3132565946033597, -0.01936272471334829, 0.13019912358366506, 0.16063533778857447, -0.011143127165783251, -0.15323703912883263, -0.3132565946033597, -0.011143127165783251, 0.4966860315069367],)

Method ArnoldiMethod.jl

Code:

using Zygote, LinearAlgebra
using ArnoldiMethod: partialschur, SR

function test_eigmin(A::AbstractMatrix)
    decomp, history = partialschur(A; nev=1, which=SR())
    ev = minimum(real.(decomp.eigenvalues))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(ArnoldiMethod._partialschur),Array{Float64,2},Type{Float64},Int64,Int64,Int64,Float64,Int64,SR}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:612
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:625
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #partialschur#1 at /Users/eh540/.julia/packages/ArnoldiMethod/5fDBS/src/run.jl:103 [inlined]
 [17] (::typeof((#partialschur#1)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #partialschur at ./none:0 [inlined]
 [19] (::typeof((#partialschur)))(::Tuple{NamedTuple{(:Q, :R, :eigenvalues),Tuple{Nothing,Nothing,Array{Float64,1}}},Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:5 [inlined]
 [21] (::typeof((test_eigmin)))(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [23] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [24] top-level scope at show.jl:582
 [25] include at ./boot.jl:328 [inlined]
 [26] include_relative(::Module, ::String) at ./loading.jl:1094
 [27] include(::Module, ::String) at ./Base.jl:31
 [28] exec_options(::Base.JLOptions) at ./client.jl:295
 [29] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/ArnoldiMethod.jl:13
eigmin(A) = -1.0756568983042305
test_eigmin(A) = -1.0756568983042305

Method GenericLinearAlgebra.jl

Code:

using Zygote, LinearAlgebra
using GenericLinearAlgebra: _eigvals!


function test_eigmin(A::AbstractMatrix)
    vals = _eigvals!(copy(A))
    ev = minimum(real.(vals))
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: TypeError: in ccall: first argument not a pointer or valid constant expression, expected Ptr, got Tuple{Symbol,String}
Stacktrace:
 [1] _forward(::Zygote.Context, ::typeof(axpy!), ::Int64, ::Float64, ::Ptr{Float64}, ::Int64, ::Ptr{Float64}, ::Int64) at /Users/julia/buildbot/worker/package_macos64/build/usr/share/julia/stdlib/v1.3/LinearAlgebra/src/blas.jl:455
 [2] _forward(::Zygote.Context, ::typeof(axpy!), ::Float64, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}, ::SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}) at ./gcutils.jl:91
 [3] lmul! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/householder.jl:64 [inlined]
 [4] _forward(::Zygote.Context, ::typeof(lmul!), ::Adjoint{Float64,GenericLinearAlgebra.Householder{Float64,SubArray{Float64,1,Array{Float64,2},Tuple{UnitRange{Int64},Int64},true}}}, ::SubArray{Float64,2,Array{Float64,2},Tuple{UnitRange{Int64},UnitRange{Int64}},false}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [5] _hessenberg! at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:53 [inlined]
 [6] _forward(::Zygote.Context, ::typeof(GenericLinearAlgebra._hessenberg!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [7] #_schur!#21 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [8] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_schur!#21")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(GenericLinearAlgebra._schur!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [9] #_eigvals!#23 at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:164 [inlined]
 [10] _forward(::Zygote.Context, ::getfield(GenericLinearAlgebra, Symbol("##_eigvals!#23")), ::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::typeof(_eigvals!), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [11] test_eigmin at /Users/eh540/.julia/packages/GenericLinearAlgebra/Lme7U/src/eigenGeneral.jl:241 [inlined]
 [12] _forward(::Zygote.Context, ::typeof(test_eigmin), ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [13] _forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:31
 [14] forward(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:37
 [15] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:46
 [16] top-level scope at show.jl:582
 [17] include at ./boot.jl:328 [inlined]
 [18] include_relative(::Module, ::String) at ./loading.jl:1094
 [19] include(::Module, ::String) at ./Base.jl:31
 [20] exec_options(::Base.JLOptions) at ./client.jl:295
 [21] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/GenericLinearAlgebra.jl:14
eigmin(A) = -1.3662825244330794
test_eigmin(A) = -1.3662825244330796

Method IterativeSolvers.jl

Code:

using Zygote, LinearAlgebra
using IterativeSolvers: lobpcg


function test_eigmin(A::AbstractMatrix)
    results = lobpcg(A, false, 1)
    ev = minimum(results.λ)
    return ev
end
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{getfield(IterativeSolvers, Symbol("##lobpcg!#53")),Bool,Int64,Bool,Float64,typeof(IterativeSolvers.lobpcg!),IterativeSolvers.LOBPCGIterator{false,Float64,Array{Float64,2},Nothing,Array{Float64,1},Array{Float64,1},Array{Int64,1},Array{Float64,2},IterativeSolvers.Blocks{false,Float64,Array{Float64,2}},IterativeSolvers.CholQR{Array{Float64,2}},IterativeSolvers.RPreconditioner{Nothing,Float64,Array{Float64,2}},IterativeSolvers.Constraint{Nothing,Array{Nothing,2},Array{Nothing,2},Nothing},IterativeSolvers.BlockGram{false,Array{Float64,2}},Array{Bool,1},Array{IterativeSolvers.LOBPCGState{Array{Float64,1},Array{Float64,1}},1}}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for at ./array.jl:612 [inlined]
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:622
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #lobpcg! at ./none:0 [inlined]
 [17] (::typeof((#lobpcg!)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] #lobpcg#52 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:837 [inlined]
 [19] (::typeof((#lobpcg#52)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #lobpcg at ./none:0 [inlined]
 [21] (::typeof((#lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] #lobpcg#50 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [23] (::typeof((#lobpcg#50)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:791 [inlined]
 [25] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] #lobpcg#49 at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [27] (::typeof((#lobpcg#49)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [28] lobpcg at /Users/eh540/.julia/packages/IterativeSolvers/QuajJ/src/lobpcg.jl:788 [inlined]
 [29] (::typeof((lobpcg)))(::NamedTuple{(:λ, :X, :tolerance, :residual_norms, :iterations, :maxiter, :converged, :trace),Tuple{Array{Float64,1},Nothing,Nothing,Nothing,Nothing,Nothing,Nothing,Nothing}}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [30] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:6 [inlined]
 [31] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [32] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [33] top-level scope at show.jl:582
 [34] include at ./boot.jl:328 [inlined]
 [35] include_relative(::Module, ::String) at ./loading.jl:1094
 [36] include(::Module, ::String) at ./Base.jl:31
 [37] exec_options(::Base.JLOptions) at ./client.jl:295
 [38] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/IterativeSolvers.jl:14
eigmin(A) = -0.8562133850903744
test_eigmin(A) = -0.8562133850628105

Method KrylovKit.jl

Code:

using Zygote, LinearAlgebra
using KrylovKit: eigsolve

test_eigmin(A::AbstractMatrix) = eigsolve(A, 1, :SR )[1][1]
A = rand(4,4); A = A + A';

@show eigmin(A)
@show test_eigmin(A)
@show Zygote.gradient(test_eigmin, A)

Results:

ERROR: LoadError: Compiling Tuple{typeof(eigsolve),Array{Float64,2},Array{Float64,1},Int64,Symbol,KrylovKit.Lanczos{KrylovKit.ModifiedGramSchmidt2,Float64}}: DimensionMismatch("dimensions must match")
Stacktrace:
 [1] promote_shape at ./indices.jl:154 [inlined]
 [2] _promote_shape at ./iterators.jl:317 [inlined]
 [3] axes at ./iterators.jl:316 [inlined]
 [4] map at ./tuple.jl:140 [inlined]
 [5] axes at ./iterators.jl:316 [inlined]
 [6] _array_for(::Type{Array{IRTools.Variable,1}}, ::Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}}, ::Base.HasShape{1}) at ./array.jl:612
 [7] collect(::Base.Generator{Base.Iterators.Zip{Tuple{Array{Any,1},Base.Iterators.Zip{Tuple{Array{Any,1},Array{Any,1}}}}},getfield(IRTools, Symbol("##126#132"))}) at ./array.jl:625
 [8] prune!(::IRTools.IR) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/passes/passes.jl:121
 [9] |> at ./operators.jl:854 [inlined]
 [10] #IR#23(::Bool, ::Bool, ::Type{IRTools.IR}, ::IRTools.Meta) at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:195
 [11] IR at /Users/eh540/.julia/packages/IRTools/y0Ot8/src/ir/wrap.jl:190 [inlined]
 [12] _lookup_grad(::Type) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/emit.jl:117
 [13] #s59#1241 at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:20 [inlined]
 [14] #s59#1241(::Any, ::Any, ::Any) at ./none:0
 [15] (::Core.GeneratedFunctionStub)(::Any, ::Vararg{Any,N} where N) at ./boot.jl:524
 [16] #eigsolve#41 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:163 [inlined]
 [17] (::typeof((#eigsolve#41)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [18] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:149 [inlined]
 [19] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [20] #eigsolve#39 at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [21] (::typeof((#eigsolve#39)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [22] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [23] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [24] eigsolve at /Users/eh540/.julia/packages/KrylovKit/kSkwM/src/eigsolve/eigsolve.jl:144 [inlined]
 [25] (::typeof((eigsolve)))(::Tuple{Array{Float64,1},Nothing,Nothing}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface2.jl:0
 [26] test_eigmin at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:4 [inlined]
 [27] (::getfield(Zygote, Symbol("##32#33")){typeof((test_eigmin))})(::Float64) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:38
 [28] gradient(::Function, ::Array{Float64,2}) at /Users/eh540/.julia/packages/Zygote/fuj2C/src/compiler/interface.jl:47
 [29] top-level scope at show.jl:582
 [30] include at ./boot.jl:328 [inlined]
 [31] include_relative(::Module, ::String) at ./loading.jl:1094
 [32] include(::Module, ::String) at ./Base.jl:31
 [33] exec_options(::Base.JLOptions) at ./client.jl:295
 [34] _start() at ./client.jl:468
in expression starting at /Users/eh540/Dropbox (Personal)/Research/Code/julia-projects/DiffEvals/Methods/KrylovKit.jl:9
eigmin(A) = -0.6743676749679843
test_eigmin(A) = -0.6743676749679844
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment