General eigen solution for sparse matrices

I have an eigen equation in the form: Ax=vBx, where A and B are large sparse matrices. In Matlab, I would use eigs(A,B,k,sigma) to solve this. k is the number of eigenvalues based on the value of sigma

How can I solve this in Julia?

Pretty much the same - see eigs with using LinearAlgebra

Many thanks. I couldnā€™t find eigs in LinearAlgebra though. eigs is not in the standard document:
https://docs.julialang.org/en/v1.0.0/stdlib/LinearAlgebra/

duplicate because my reply was stuck in new user limbo

Looks like it is in the Arpack package now.

I ended up installing julia 0.7 in parallel with julia 1.0 because of all these changes since 0.6. At least 0.7 tells you where those functions went or what function replaces them.

Coming from matlab, julia 0.6 was rather easy to transition to. So many functions in the Base package with identical names and very close or identical syntax. It looks like it is not the case any more in 0.7 / 1.0. Some basic functions (for a matlab user) have moved to dedicated packages and some others have been replaced by other functions / been renamed. It may be for the best in the long term but it will certainly make it more difficult for matlab users to transition to julia in the short term. Unless they use julia 0.7 for the transition. Maybe.

Exactly. I have been using Matlab extensively, and Matlab has well-established user manual and documentation. For Julia, it is not easy to search for something that I am trying to use. The standard library introduces basic eigen equation solutions, but not for large and sparse matrices. If it is in the Arpack package, it shouldnā€™t be difficult to incorporate.

I had a go with Arpack.jl, and the file is not working for me:
https://github.com/JuliaLinearAlgebra/Arpack.jl/commit/9da9aa0a0aaafdcca205e0829418a09f08d72d10

It seems that I am not the only person that has this problem.

The julia documentation and packages certainly need some work but this is going according to plan apparently. Julia 1.0 only signals that the language will not break compatibility in the future. It does not mean that the documentation and standard packages are ready for professional use.

This is unfortunate as many users may want to try Julia now thanks to the landmark 1.0 announcement, only to be disappointed by the current state of the ecosystem. A failed opportunity to kickstart Julia adoption in my opinion.

A better signal for when Julia will be ready for production may be when JuliaPro transitions to 1.x (or maybe 0.7.x). I doubt this annoucement will be seen by as many people as the one for Julia 1.0 though.

1 Like

If you want, you can try out IterativeSolvers.jl which has the LOBPCG algorithm:

julia> using IterativeSolvers

julia> using SparseArrays, LinearAlgebra

julia> A = sprand(100,100,0.1); A = A + A' + 10I;

julia> B = sprand(100,100,0.1); B = B + B' + 20I;

julia> largest = true;

julia> nev = 1;

julia> r = lobpcg(A, B, largest, nev)
Results of LOBPCG Algorithm
 * Algorithm: LOBPCG - CholQR
 * Ī»: [0.8138294334999978]
 * Residual norm(s): [1.794873755112849e-5]
 * Convergence
   * Iterations: 29
   * Converged: true
   * Iterations limit: 200

julia> r.X, r.Ī»
([-0.0222962; -0.00650525; ā€¦ ; 0.00265957; 0.0222808], [0.813829])

I am assuming your matrices are Hermitian and B is positive definite.

1 Like

It is a little disappointing that some basic functions are not in place yet. However, I know that it is a relatively new language, and it would take some time to get to know a new language. I am just trying to see if it is possible to move from Matlab to Julia at the moment, and it seems that it is not so easy.

Many thanks. It looks like a nice function to try out. Although I couldnā€™t add this package yet. I tried ā€˜Pkg.add(ā€œIterativeSolversā€)ā€™, and I got an ā€˜ERROR: UndefVarError: Pkg not definedā€™. Not sure what is wrong yet.

Sometimes my matrices are Hermitian matrices (positive definite), but not always. It takes more effort to construct equations so that matrices are Hermitian, and sometimes not feasible.

You need to first import Pkg to load the package manager.

Ok. Many thanks. IterativeSolvers.jl is working for me now. A nice alternative to eigs !

I have added ā€˜Arpackā€™ as well. Maybe ā€˜eigsā€™ is also working in Julia. I need to test these functions.

On the contrary, switching from Matlab to Julia is relatively easy. It is true that finding out the mapping from Matlab functionality to Julia packages takes a bit of work, but the Discourse community is an amazing resource: just ask a question and you will get it answered typically within minutes. To port individual functions is made relatively straightforward by the one-based indexing in both languages.

For FEA functionality that you were asking about above you could also check out https://github.com/PetrKryslUCSD/FinEtools.jl: For instance the generalized eigenvalue problem is solved in many examples (vibration).

4 Likes

Ok. Sounds very encouraging. The reason that I raised this question is that eigen solutions can be used to calculate natural frequencies and mode shapes etc, These are basic finite element models, yet widely used in some applications. For some problems, A and B matrices are both real and symmetric. However, for some complicated problems, A and B are complex and non-Hermitian, which makes a ā€˜eigsā€™ like function necessary.

Thanks for the FinEtools link. I will look into this.

May I ask what kind of applications are these?

These are from guided wave applications, where the wave guide is considered to be infinitely long and analytical solutions have to be introduced in the wave propagation direction. This would introduce complex solutions into the equations. Non-Hermitian structure could come from coupling between different regions in the waveguide.

2 Likes