Zero pivot in EUCLID preconditioner

[valonbb]

I am aware that EUCLID is no longer supported by hypre team, but I find it useful to use it as a preconditioner and in particular with IJ interface. I have realised that EUCLID does not support matrices with zero pivoting, so is there any way around it to activate that or what would it be the best approach around it? Thank you.

[oseikuffuor1]

@valonbb have you tried the HYPRE_ILU options? This is the replacement for EUCLID and I believe some of the options have support for zero diagonals. Alternatively, you can perturb your input matrix with small entries corresponding to the zero pivots, prior to calling EUCLID.

[valonbb]

@oseikuffuor1 Initially that is what I have done, used HYPRE_ILU, rather than EUCLID, but I have found its behaviour to not be very robust with varying number of mpi processes. Replacing zero pivots with small entries sound like a good idea, what range of "small numbers" (relative to drop tol or min/max of matrix at hand) would you recommend to use to replace zero pivots?

Thanks a lot for all your help.

[oseikuffuor1]

@valonbb something small, relative to the entries in the matrix should be good. droptol * min/max is fine if you know what min and max are. You can also pick something a few orders of magnitude smaller than your min (nonzero) diagonal entry.

Have you tried any reordering techniques (RCM, AMD, etc)? If you have access to the matrix info, you could also keep the zero (or small) diagonals close to the lower part of the matrix to increase the chances of nonzero pivots during the factorization.

[valonbb]

@oseikuffuor1 Thank you for your suggestions. I have not actually tried using any reordering techniques, as it does not seem that EUCLID supports that explicitly, but sure that can be done before EUCLID is called.

you could also keep the zero (or small) diagonals close to the lower part of the matrix to increase the chances of nonzero pivots during the factorization.

Yes, I have access to matrix. Would this mean some further permutation of matrix to push smaller diagonal values on the lower part? I am not sure if I understood this correctly.

Thanks a lot for all your help, very much appreciated

[oseikuffuor1]

@valonbb Yes, if you have access to the matrix, then permuting like you said could help. Keep in mind that depending on the structure of the matrix, that will also affect the fill-in pattern of the resulting factorization. I assume if you are using EUCLID then you are relying on ILU(k) routines? Also, did perturbing small entries not help?

[valonbb]

@oseikuffuor1 Yes perturbing small entries actually help and it does seem to work well so far - thanks a lot for your suggestion and all your help. You are certainly right, this would depend on the structure of the matrix and so affecting fill-ins in the factorization. Yes, I am relying on ILUK as I think EUCLID supports ILUT in the sequential operation, only?

[oseikuffuor1]

@valonbb Yes EUCLID's ILUT can only be used sequentially, unfortunately. Curious, are you using the block Jacobi option for EUCLID or the default parallel ILU option? I am still surprised that HYPRE_ILU does not work for your problem. Are you able to share what the structure of your system looks like?

[valonbb]

Hi @oseikuffuor1,

Thanks a lot for your reply and apologise that it is taking this long from my side to reply.

I am mainly using EUCLID default option PILU, rather than block Jacobi. Although, the pivot is now fixed based on your suggestions – thanks a lot for your help, I find EUCLID to be relatively slow compare to ILU (sequentially).

I am indeed able to share matrix structure with you, please find attached the matrix in IJ format and matrix market format, hopefully this helps.
A.txt - matrix market format
IJ.A.txt - IJ format

For this problem I have used BiCGSTAB solver and ILU preconditioner along with these parameters:

HYPRE_BiCGSTABSetMaxIter(solver, 1000);
HYPRE_BiCGSTABSetTol(solver, 1,0e-10);

HYPRE_ILUSetMaxIter(precond, 1);
HYPRE_ILUSetTol(precond, 0.);
HYPRE_ILUSetType(precond, 31); /* I have also used other options here: 0, 1, 10, 11, 31 */
HYPRE_ILUSetDropThreshold(precond, 1.0e-3);

When I run this sequentially, that is using mpirun -np 1 – I get
number of iterations = 9
final residual norm = 1.98434e-11

Whereas using 2 processes – mpirun -np 2
number of iterations: 1000
final residual norm: 0.99994
Suggesting that solution has not converged. I find similar behaviour with changing the type of ilu and also varying number of processors. Also, increasing maximum number of iterations does not seem to help.

I very much appreciate your time to look after this - thanks a lot for all your help.

[oseikuffuor1]

@valonbb thanks for the data. Can you also share the data for the two processor case? You can use hypre's IJMatrixPrint routine and it should print the matrix decomposed as you have it on the two processors (so two systems, one for each processor).

[valonbb]

@oseikuffuor1 yes, please find below the matrix decomposed across processors:
IJ.out.A.00001.txt
IJ.out.A.00000.txt

Thanks a lot for your help.