wilkoj
wilkoj

Reputation: 23

SIGSEV fault in scipy's eigh

I am having problems with a 'SIGSEGV' fault when diagonalising matrices with scipy's linalg.eigh.

The following code reproduces the problem after a few iterations of the for loop:

import scipy.sparse as sparse
import scipy.linalg as linalg

a = sparse.random(128, 128, 0.8)
b = sparse.random(128, 128, 0.8)*1j
a = a + b
a = a + a.H

for i in range(0, 100):
        E, R = linalg.eigh(a.todense())

After stepping and debugging, the faults seem to happen when scipy calls LAPACK, but I am unsure as to what I can do to stop this from happening.

I am running python 3.7, using scipy 1.2.1 under macOS Catalina.

Any help would be much appreciated.

Upvotes: 2

Views: 95

Answers (1)

v-joe
v-joe

Reputation: 366

The problem is related to not all memory being freed that is used for E and R and new E and R arrays being created every iteration. Also inserting a gc.collect() at end of the loop doesn't fix this. Deleting E and R after the end of each iteration, or making sure to not creating new array copies for E and R in every iteration fixes the issue, as shown in the two code snippets:

#manually delete E and R at the end of each iteration
import scipy.sparse as sparse
import scipy.linalg as linalg

a = sparse.random(128, 128, 0.8)
b = sparse.random(128, 128, 0.8)*1j
a = a + b
a = a + a.H

for i in range(0, 100):
        E, R = linalg.eigh(a.todense())
        del(E)
        del(R)
#create E and R before the loop and store eigenpairs directly into these arrays
import scipy.sparse as sparse
import scipy.linalg as linalg
import numpy as np

a = sparse.random(128, 128, 0.8)
b = sparse.random(128, 128, 0.8)*1j
a = a + b
a = a + a.H

E = np.empty(128, np.float)
R = np.empty((128, 128), np.complex)

for i in range(0, 100):
        E[:], R[:][:] = linalg.eigh(a.todense())

Upvotes: 3

Related Questions