Reduce memory usage for locking==false

[twesterhout]

In Section 3.4.1 (Computational Requirements: Memory) of your paper "PRIMME: PReconditioned Iterative MultiMethod Eigensolver: Methods and software description" it is mentioned that the amount of allocated memory can be reduced when locking is disabled.

In the code I see that V and W are allocated unconditionally, so it seems like there is currently no way to reduce the amount of used memory.

When diagonalizing very big matrices (which are then never actually stored in memory), getting rid of numEvals vectors is a very desirable feature. Are there any plans to implement it? Would you accept PRs in this direction?

Thanks!

[eromero-vlc]

I agree that saving space is a desirable feature, and we didn't put the feature in PRIMME 3 because the memory savings are moderate (less than a factor of two). Nevertheless PR are welcome and I am willing to put effort in integrating them into the project!

One easy way that I see of how to do this is by adding a flag in primme_params indicating whether the user has passed a big enough space in evecs to hold maxBasisSize vectors. If that's the case, main_iter will use as V the pointer provided as evecs instead of allocating a new space.

About W, PRIMME stores the A*V vectors in W, where A is the matrix to diagonalize. Unfortunately saving the space of W is more challenging and will require extra computations and a deeper knowledge of how eigensolvers and PRIMME work.

[twesterhout]

About W, PRIMME stores the A*V vectors in W, where A is the matrix to diagonalize. Unfortunately saving the space of W is more challenging and will require extra computations and a deeper knowledge of how eigensolvers and PRIMME work.

From Algorithm 3.1 I get the impression that the size of W can be reduced from (n, m_max) to (n, b) (i.e. from maxBasisSize to maxBlockSize columns). H = V^T * W calculation has to be split into m_max / b steps (skipping cache locality implications of it for now), but no additional computations are necessary? or am I missing something?

[eromero-vlc]

Notice that each loop of the while starting at step 5 updates b columns of V (let's ignore step 6). Updating H (step 13) can be done with W storing the product of A times the last b columns of V. But computing w^(m) = W * s_{0:b-1} (step 7), needs the whole W.

It is possible to do step 7 like this w^(m) = A * (V * s_{0:b-1}), which is, first, compute V * s_{0:b-1}, then multiply it by A, and finally store the result on w^(m). But this approach require to multiple by A more times and that can be more expensive sometimes.