I have a sparse, structurally-symmetric matrix. The solution of this matrix is part of a larger iterative solver, so the matrix needs to be solved many times. The values in the matrix will modestly change between iterations, but the sparsity pattern will remain the same.
Based on my reading of the PARDISO documentation, it sounds like this situation argues strongly for the implementation of the CGS iterative scheme. It sounds like PARDISO needs to be called once with phase=11, and then subsequently called with phase=23 and iparm(4)=31. Is this correct? I havent been able to get this to work. I get the following errors:
*** Error in PARDISO ( numerical_factorization) error_num= -3
Are there any public examples that would be useful? Any other thoughts?
Thanks for looking into this. You are correct that I was previously using mtype=1. Using mtype=11 appears to produce the same errors.
A test case is attached. I suspect that there may be something basic wrong with my input. Please let me know what you find.
I've re-attached the file to this posting, updated with the values of the RHS that I'm using. (The same RHS is used for each matrix.)
However, I suspect the problem I'm having is related to how I'm calling the subroutine. I have no problem using the "default" iparm values to call PARDISO and obtain a solution to the matrices. That works just fine. When I try to use the CGS iteration (iparm(4)=31 in Fortran, iparm=31in C) to replace the computation of LU is when I encounter problems.
Based on my reading of the documentation, my best guess as to how this is supposed to work is as follows (example in abbreviated Fortran):
sparse_values = ...
iparm(4) = 0
When I attempt this, I get the following:
At what point in the pardiso_unsym_c are you modifying iparm(3)?
After eliminating every other possibility,I noticed that I was using single-precision arrays whereas your code was using double-precision arrays. When I switched my arrays to double-preicision, the problem went away.
Since the matrices aren't going to be horribly large, I suppose we can use double-preicision arrays, but is there some reason single-precision isn't working? Could it be made to work?
If Line 6 is commented out, Line 7 in un-commented, and Line 50 is commented out, the program will execute correctly. I'm using 12.0.3.