# What range of LinearAlgebra does Enzyme support?

**URL:** <https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552>\
**Category:** Numerics\
**Tags:** question, differentiation, enzyme\
**Created:** [March 13, 2024, 6:18am UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552 "2024-03-13T06:18:11Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![ToPo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/topo/32/205868_2.png) [@ToPo](https://discourse.julialang.org/u/ToPo)\
**Post date:** [March 13, 2024, 6:18am UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/1 "2024-03-13T06:18:11Z")

</div>

Hi all,

I would like to use the LinearAlgebra package in Enzyme.  
Since Julia’s sparse matrix package, SparseArrays, did not seem to be available, I created my own sparse-matrix library.

When I looked at the Enzyme testing program, I found a BLAS test.  
Therefore, using that library, I created a program for the ICCG method to solve the simultaneous equations.  
I used BLAS.dot(), but in some cases the program crashed when used.  
In the program, resl0 = sqrt(BLAS.dot(res, res)) is not a problem. But sk = BLAS.dot(res, zk2) and pkapk = BLAS.dot(pk, zk) crash the program.

I would like to know if Enzyme currently supports up to the LinearAlgebra part.  
I am not even able to use basic functions such as det(), although I have not used them in this ICCG program.  
If you know more about this, please let me know.

Some of the programs are as follows:  
#------------------------------------------------------------------------------------

# 

#------------------------------------------------------------------------------------  
function ICCG(A::MySparseMatrixCRS, b::Vector{Float64}, x0::Vector{Float64}=zeros(Float64, length(b)), itr\_max::Int64=100000, tol\_abs::Float64=1.0e-13, tol\_rel::Float64=1.0e-13)  
#  
small = 1.0e-15

```
# Initalize working value and arrays
n = A.nrows
pk = zeros(Float64, n)
x = zeros(Float64, n)

# Calculate initial residual vector and the norm
res = copy(b) - A * x0

# L^1-norm of residual
#resl0 = sqrt(my_dot_product(res, res))
resl0 = sqrt(BLAS.dot(res, res))
tol = tol_abs + tol_rel * resl0

# Initially converged solution
if resl0 < tol
    println("Iteration: 0, Final residual = $(resl0)")
    return x
end

# create diag
diag = Vector{Int64}(undef, n)
Threads.@threads for i in 1 : n
     for k in A.row_ptr[i] : A.row_ptr[i+1] - 1
         if A.col_ind[k] == i
             diag[i] = k
        end
    end
end

# Calculate elements of diagonal preconditioning matrix
d = Vector{Float64}(undef, n)
for i in 1 : n
     d[i] = A.values[diag[i]]
     for k in A.row_ptr[i] : diag[i] - 1
         d[i] -= A.values[k]^2.0 * d[A.col_ind[k]]
    end
     d[i] = 1.0 / (d[i] + small)
end

# Start iterations
s0 = 1.0e+20
itr_used = 0
 for l in 1 : itr_max
    # forward substitution
    zk1 = forward_substitution(A, d, diag, res)

    # Backward substitution
    zk2 = backward_substitution(A, d, diag, zk1)

    # Inner product
    sk = my_dot_product(res, zk2)
    #sk = BLAS.dot(res, zk2)

    # Calculate beta
    bet = sk / s0

    # Calculate new search vector pk
    pk .= zk2 .+ bet .* pk

    # Calculate scalar product (pk.a pk) and alpha (overwrite zk)
    zk = A * pk

    # Inner product
    pkapk = my_dot_product(pk, zk)
    #pkapk = BLAS.dot(pk, zk)

    # Calculate alpha
    alf = sk / pkapk

    # Update solution vector
    x .+= alf .* pk

    # Update residual vector
    res .-= alf .* zk

    # L1-norm of residual
    #resl = my_dot_product(res, res)
    resl = BLAS.dot(res, res)

    # update
    s0 = sk
    itr_used = itr_used + 1

    # Check convergence
    rsm = sqrt(resl)
    if rsm < tol
        #println("Convergence achieved after $(itr_used) iterations, residual = $(rsm).")
        return x
    end
end

#
println("Max iteration. Not converged")
return x

```

end

---

<div class="post-metadata">

**Author:** ![wsmoses](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/wsmoses/32/26497_2.png) [@wsmoses](https://discourse.julialang.org/u/wsmoses)\
**Post date:** [March 13, 2024, 7:34am UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/2 "2024-03-13T07:34:57Z")

</div>

gemm (matmul), gemv(matvec), dot (dot product), axpy (vector combine), and others exist with fast blas libraries for float32/64. Fallback derivatives for other BLAS libraries is also included. Lapack support (with the exception of some additions and current work underway, such as cholesky, etc) is presently limited. I believe det is a lapack call in julia. It shouldn’t be hard to add if you’re interested.

As for whatever else you saw, it would be hard to help without seeing the error message / backtrace.

---

<div class="post-metadata">

**Author:** ![ToPo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/topo/32/205868_2.png) [@ToPo](https://discourse.julialang.org/u/ToPo)\
**Post date:** [March 13, 2024, 2:03pm UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/3 "2024-03-13T14:03:23Z")

</div>

For example, the following changes were made in the ICCC program:  
sk = my\_dot\_product(res, zk2) → sk = BLAS.dot(res, zk2)

my\_dot\_product() is a function of the linear algebra library I created in Julia.  
The types of “res” and “zk2” are both Vector{Float64}.

When run with the normal REPL, Julia itself crashes and stops before checking for errors.  
Therefore, we ran it in debug mode.  
Sorry it is very long. Here are the errors.

* * *

[51880] signal (11.2): Segmentation fault: 11  
in expression starting at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/debugger/run\_debugger.jl:13  
\_ZN13GradientUtils15cacheForReverseERN4llvm9IRBuilderINS0\_14ConstantFolderENS0\_24IRBuilderDefaultInserterEEEPNS0\_5ValueEib at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator10handle\_dotE8BlasInfoRN4llvm8CallInstEPNS1\_8FunctionERKNSt3\_\_16vectorIbNS6\_9allocatorIbEEEEPNS1\_4TypeE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator10handleBLASERN4llvm8CallInstEPNS0\_8FunctionE8BlasInfoRKNSt3\_\_16vectorIbNS6\_9allocatorIbEEEE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
_ZN16AdjointGenerator26handleKnownCallDerivativesERN4llvm8CallInstEPNS0\_8FunctionENS0\_9StringRefERKNSt3\_\_16vectorIbNS6\_9allocatorIbEEEEPS1_ at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator13visitCallInstERN4llvm8CallInstE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN11EnzymeLogic23CreatePrimalAndGradientE14RequestContextOK15ReverseCacheKeyR12TypeAnalysisPK15AugmentedReturnb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator28recursivelyHandleSubfunctionERN4llvm8CallInstEPNS0\_8FunctionERKNSt3\_\_16vectorIbNS5\_9allocatorIbEEEEb10DIFFE\_TYPEb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator13visitCallInstERN4llvm8CallInstE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN11EnzymeLogic23CreatePrimalAndGradientE14RequestContextOK15ReverseCacheKeyR12TypeAnalysisPK15AugmentedReturnb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator28recursivelyHandleSubfunctionERN4llvm8CallInstEPNS0\_8FunctionERKNSt3\_\_16vectorIbNS5\_9allocatorIbEEEEb10DIFFE\_TYPEb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator13visitCallInstERN4llvm8CallInstE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN11EnzymeLogic23CreatePrimalAndGradientE14RequestContextOK15ReverseCacheKeyR12TypeAnalysisPK15AugmentedReturnb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator28recursivelyHandleSubfunctionERN4llvm8CallInstEPNS0\_8FunctionERKNSt3\_\_16vectorIbNS5\_9allocatorIbEEEEb10DIFFE\_TYPEb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN16AdjointGenerator13visitCallInstERN4llvm8CallInstE at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
\_ZN11EnzymeLogic23CreatePrimalAndGradientE14RequestContextOK15ReverseCacheKeyR12TypeAnalysisPK15AugmentedReturnb at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
EnzymeCreatePrimalAndGradient at /Users/iberico\_mac/.julia/artifacts/f9ff3c847f566c574f14defcc79f62b90d0c2cc8/lib/libEnzyme-15.dylib (unknown line)  
EnzymeCreatePrimalAndGradient at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/api.jl:154  
unknown function (ip: 0x28c54c25b)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
enzyme! at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:3113  
unknown function (ip: 0x28c4fb153)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
#codegen#488 at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:4825  
codegen at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:4340 [inlined]  
\_thunk at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:5430  
\_thunk at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:5430 [inlined]  
cached\_compilation at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:5464 [inlined]  
#532 at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:5530  
JuliaContext at /Users/iberico\_mac/.julia/packages/GPUCompiler/U36Ed/src/driver.jl:47  
unknown function (ip: 0x281030263)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
#s1883#531 at /Users/iberico\_mac/.julia/packages/Enzyme/wR2t7/src/compiler.jl:5482 [inlined]  
#s1883#531 at ./none:0  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
GeneratedFunctionStub at ./boot.jl:602  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
jl\_call\_staged at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/method.c:540  
ijl\_code\_for\_staged at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/method.c:593  
get\_staged at ./compiler/utilities.jl:123  
get\_staged at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/utils.jl:216 [inlined]  
#prepare\_framecode#40 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/construct.jl:157  
prepare\_framecode at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/construct.jl:138 [inlined]  
#prepare\_call#43 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/construct.jl:249  
prepare\_call at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/construct.jl:235 [inlined]  
#get\_call\_framecode#58 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/localmethtable.jl:63  
get\_call\_framecode at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/localmethtable.jl:10 [inlined]  
#evaluate\_call\_recurse!#65 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:274  
evaluate\_call\_recurse! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:249 [inlined]  
eval\_rhs at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:422  
unknown function (ip: 0x12549e8bb)  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:495  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:633  
finish! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:14  
finish\_and\_return! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:30  
#evaluate\_call\_recurse!#65 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:288  
evaluate\_call\_recurse! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:249 [inlined]  
eval\_rhs at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:422  
unknown function (ip: 0x12549e8bb)  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:578  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:633  
finish! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:14  
finish\_and\_return! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:30  
#evaluate\_call\_recurse!#65 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:288  
evaluate\_call\_recurse! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:249 [inlined]  
eval\_rhs at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:422  
unknown function (ip: 0x12549e8bb)  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:578  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:633  
finish! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:14  
finish\_and\_return! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:30  
#evaluate\_call\_recurse!#65 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:288  
evaluate\_call\_recurse! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:249 [inlined]  
eval\_rhs at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:422  
unknown function (ip: 0x12549e8bb)  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:578  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:633  
finish! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:14  
finish\_and\_return! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:30  
#evaluate\_call\_recurse!#65 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:288  
evaluate\_call\_recurse! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:249 [inlined]  
eval\_rhs at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:422  
unknown function (ip: 0x12549e8bb)  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:573  
step\_expr! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/interpret.jl:633  
finish! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:14  
finish\_and\_return! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:30  
unknown function (ip: 0x1254c409b)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
finish\_stack! at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:60  
unknown function (ip: 0x1255d4a5b)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
#debug\_command#85 at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:516  
debug\_command at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/JuliaInterpreter/src/commands.jl:448  
unknown function (ip: 0x1255cc06f)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
jl\_apply at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/./julia.h:1982 [inlined]  
jl\_f\_\_call\_latest at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/builtins.c:812  
#invokelatest#2 at ./essentials.jl:892 [inlined]  
invokelatest at ./essentials.jl:889 [inlined]  
our\_debug\_command at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/DebugAdapter/src/debugger\_core.jl:67  
startdebug at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/DebugAdapter/src/packagedef.jl:106  
startdebugger at /Users/iberico\_mac/.vscode/extensions/julialang.language-julia-1.75.2/scripts/packages/VSCodeDebugger/src/VSCodeDebugger.jl:43  
unknown function (ip: 0x10c400563)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
jl\_apply at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/./julia.h:1982 [inlined]  
do\_call at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/interpreter.c:126  
eval\_body at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/interpreter.c:0  
jl\_interpret\_toplevel\_thunk at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/interpreter.c:775  
jl\_toplevel\_eval\_flex at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/toplevel.c:934  
jl\_toplevel\_eval\_flex at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/toplevel.c:877  
ijl\_toplevel\_eval at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/toplevel.c:943 [inlined]  
ijl\_toplevel\_eval\_in at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/toplevel.c:985  
eval at ./boot.jl:385 [inlined]  
include\_string at ./loading.jl:2076  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
\_include at ./loading.jl:2136  
include at ./Base.jl:495  
jfptr\_include\_46545 at /Users/iberico\_mac/.julia/juliaup/julia-1.10.2+0.aarch64.apple.darwin14/lib/julia/sys.dylib (unknown line)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
exec\_options at ./client.jl:318  
\_start at ./client.jl:552  
jfptr\_\_start\_82857 at /Users/iberico\_mac/.julia/juliaup/julia-1.10.2+0.aarch64.apple.darwin14/lib/julia/sys.dylib (unknown line)  
\_jl\_invoke at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:0 [inlined]  
ijl\_apply\_generic at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/gf.c:3076  
jl\_apply at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/./julia.h:1982 [inlined]  
true\_main at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/jlapi.c:582  
jl\_repl\_entrypoint at /Users/julia/.julia/scratchspaces/a66863c6-20e8-4ff4-8a62-49f30b1f605e/agent-cache/default-honeycrisp-HL2F7YQ3XH.0/build/default-honeycrisp-HL2F7YQ3XH-0/julialang/julia-release-1-dot-10/src/jlapi.c:731  
Allocations: 282966118 (Pool: 282828385; Big: 137733); GC: 205

---

<div class="post-metadata">

**Author:** ![ToPo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/topo/32/205868_2.png) [@ToPo](https://discourse.julialang.org/u/ToPo)\
**Post date:** [March 13, 2024, 2:23pm UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/4 "2024-03-13T14:23:36Z")

</div>

I have incorporated the entire finite element method program into Enzyme’s automatic differentiation.  
As a linear algebra library I need the following:  
dot(), inv(), det(), transpose() (or I want to calculate B’ \* C \* B )

The sparse arrays sparse matrix creation and direct method linear solver are perfect.  
However, it still takes time?

---

<div class="post-metadata">

**Author:** ![vchuravy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vchuravy/32/8_2.png) [@vchuravy](https://discourse.julialang.org/u/vchuravy)\
**Post date:** [March 13, 2024, 2:53pm UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/5 "2024-03-13T14:53:50Z")

</div>

Please format your messages (see [Please read: make it easier to help you](https://discourse.julialang.org/t/please-read-make-it-easier-to-help-you/14757))

Also it helps tremendously if you post minimal working examples so that we can more easier reproduce what you are seeing

---

<div class="post-metadata">

**Author:** ![ToPo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/topo/32/205868_2.png) [@ToPo](https://discourse.julialang.org/u/ToPo)\
**Post date:** [March 14, 2024, 12:36am UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/6 "2024-03-14T00:36:57Z")

</div>

Sorry about that. I will be careful.  
Also, thanks for the information on how to write.

---

<div class="post-metadata">

**Author:** ![ToPo](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/topo/32/205868_2.png) [@ToPo](https://discourse.julialang.org/u/ToPo)\
**Post date:** [March 14, 2024, 3:08am UTC](https://discourse.julialang.org/t/what-range-of-linearalgebra-does-enzyme-support/111552/7 "2024-03-14T03:08:21Z")

</div>

Once again, I have created a program to test independently.

Using dot() in the calculation of variable sk causes the program to crash.  
On the other hand, the calculation can be done well using my\_dot\_product(), a function of my own creation.

resl=dot(res, res) can be computed successfully.  
Why is it sometimes calculated well and sometimes not?

```julia
using LinearAlgebra
using Enzyme
import Base.*
#----------------------------------------------------------------------------
# sparse matrix (CRS format)
#----------------------------------------------------------------------------
mutable struct MySparseMatrixCRS
    values::Vector{Float64} # non zero value
    col_ind::Vector{Int64} # index
    row_ptr::Vector{Int64} # pointer
    nrows::Int64
    ncols::Int64
end
#----------------------------------------------------------------------------
# A * b : for my sparse matrix type
#----------------------------------------------------------------------------
function *(A::MySparseMatrixCRS, x::Vector{Float64})
    z = Vector{Float64}(undef, A.nrows)
    for i in 1 : A.nrows
        z[i] = 0.0
        for j in A.row_ptr[i] : A.row_ptr[i+1] - 1
            z[i] += A.values[j] * x[A.col_ind[j]]
        end
    end
    return z
end
#----------------------------------------------------------------------------
# dot product by myself
#----------------------------------------------------------------------------
function my_dot_product(a::Vector{Float64}, b::Vector{Float64})
    value = 0.0
    for i in 1 : length(a)
        value += a[i] * b[i]
    end
    return value
end
#----------------------------------------------------------------------------
# dense matrix -> sparse matrix
#----------------------------------------------------------------------------
function convert_to_sparse_matrix_crs(A::Matrix{Float64})
    nrows, ncols = size(A)
    values = Vector{Float64}()
    col_ind = Vector{Int64}()
    row_ptr = Vector{Int64}()
    push!(row_ptr, 1)
     for i in 1:nrows
        for j in 1:ncols
            if A[i, j] != 0.0
                push!(values, A[i, j])
                push!(col_ind, j)
            end
        end
        push!(row_ptr, length(values) + 1)
    end
    return MySparseMatrixCRS(values, col_ind, row_ptr, nrows, ncols)
end
#------------------------------------------------------------------------------------
# linear solver
#------------------------------------------------------------------------------------
function ICCG(A::MySparseMatrixCRS, b::Vector{Float64}, itr_max::Int64=100000, tol_abs::Float64=1.0e-13, tol_rel::Float64=1.0e-13)    
    #
    small = 1.0e-15
    
    # Initalize working value and arrays
    n = A.nrows
    pk = zeros(Float64, n)
    x = zeros(Float64, n)

    # Calculate initial residual vector and the norm
    res = copy(b) - A * x0

    # L^1-norm of residual
    resl0 = sqrt(dot(res, res))
    tol = tol_abs + tol_rel * resl0

    # Initially converged solution
    if resl0 < tol
        println("Iteration: 0, Final residual = $(resl0)")
        return x
    end

    # create diag
    diag = Vector{Int64}(undef, n)
    Threads.@threads for i in 1 : n
         for k in A.row_ptr[i] : A.row_ptr[i+1] - 1
             if A.col_ind[k] == i
                 diag[i] = k
            end
        end
    end

    # Calculate elements of diagonal preconditioning matrix
    d = Vector{Float64}(undef, n)
    for i in 1 : n
         d[i] = A.values[diag[i]]
         for k in A.row_ptr[i] : diag[i] - 1
             d[i] -= A.values[k]^2.0 * d[A.col_ind[k]]
        end
         d[i] = 1.0 / (d[i] + small)
    end

    # Start iterations
    s0 = 1.0e+10
    itr_used = 0
     for l in 1 : itr_max
        # forward substitution
        zk1 = forward_substitution(A, d, diag, res)

        # Backward substitution
        zk2 = backward_substitution(A, d, diag, zk1)

        # Inner product
        #sk = my_dot_product(res, zk2)
        sk = dot(res, zk2)

        # Calculate beta
        bet = sk / s0

        # Calculate new search vector pk
        pk .= zk2 .+ bet .* pk

        # Calculate scalar product (pk.a pk) and alpha (overwrite zk)
        zk = A * pk

        # Inner product
        pkapk = my_dot_product(pk, zk)
        #pkapk = dot(pk, zk)

        # Calculate alpha
        alf = sk / pkapk

        # Update solution vector
        x .+= alf .* pk

        # Update residual vector
        res .-= alf .* zk

        # L1-norm of residual
        resl = dot(res, res)

        # update
        s0 = sk
        itr_used = itr_used + 1

        # Check convergence
        rsm = sqrt(resl)
        if rsm < tol
            #println("Convergence achieved after $(itr_used) iterations, residual = $(rsm).")
            return x
        end
    end
    #
    println("Max iteration. Not converged")
    return x
end
#------------------------------------------------------------------------------------
#
#------------------------------------------------------------------------------------
function forward_substitution(A::MySparseMatrixCRS, d::Vector{Float64}, diag::Vector{Int64}, res::Vector{Float64})
    small = 1.0e-15
    zk = similar(res)
    for i in 1 : A.nrows
        zk[i] = res[i]
        for k in A.row_ptr[i] : diag[i] - 1
            zk[i] = zk[i] - A.values[k] * zk[A.col_ind[k]]
        end
        zk[i] = zk[i] * d[i]
    end
    for i in 1 : A.nrows
         zk[i] /= (d[i] + small)
    end
    return zk
end
#------------------------------------------------------------------------------------
#
#------------------------------------------------------------------------------------
function backward_substitution(A::MySparseMatrixCRS, d::Vector{Float64}, diag::Vector{Int64}, zk1::Vector{Float64})
    zk = copy(zk1)
    for i in A.nrows : -1 : 1
        for k in diag[i] + 1 : A.row_ptr[i+1] - 1
            zk[i] -= A.values[k] * zk[A.col_ind[k]]
        end
        zk[i] *= d[i]
    end
    return zk
end
#----------------------------------------------------------------------------
# test function
#----------------------------------------------------------------------------
function iccg_test(x)

    # make sparse matrix
    A_dense = Matrix{Float64}(undef, 3, 3)
    for i in 1 : 3
        for j in 1 : 3
            A_dense[i, j] = 2.0 * i + 3.0 * j^2.0 + x[i]^2.0 + x[j]^3.0
        end
    end
    # for symmetric matrix
    for i in 1 : 3
        for j in i : 3
            A_dense[i, j] = A_dense[j, i]
        end
    end
    A_sparse = convert_to_sparse_matrix_crs(A_dense)

    # right hand side
    b = Vector{Float64}(undef, 3)
    b[1] = 1.0
    b[2] = 0.0 
    b[3] = x[2]^3.0

    # solve
    u = ICCG(A_sparse, b)

    # for scalar
    return dot(u, u)
end
#----------------------------------------------------------------------------
# main
#----------------------------------------------------------------------------
x0 = rand(Float64, 3)

val = iccg_test(x0)
display(val)

df = gradient(Reverse, iccg_test, x0)
display(df)

```