From: Rick F. <ri...@ar...> - 2020-02-04 23:00:36
|
Hello all, I am solving large sparse symmetric systems. As a result I have been using the CG solver. However the only preconditioner that seems to work robustly across different matrices is the row_scaling preconditioner. I have tried jacobi_precond and chow_patel_icc_precond with the CG solver and while they work for some symmetric problems, they either fail to converge or give NaNs for other symmetric problems. The same problem converges just fine without the preconditioners or when using the row_scaling preconditioner. Are there additional special settings beyond what is provided in the documentation needed to successfully use these preconditioners? Are there other limitations of these preconditioners? My matrices regularly contain non-zero entries O(10^6). I have verified they are symmetric. I am using an NVIDIA GeForce GTX 1080 Ti GPU and ViennaCL 1.7.1. Exemplary code is: vector<map<unsigned int, double>> AA; /* Sparse matrix. */ vector<double> b; /* Right hand side vector. */ // Assemble AA and b here // ... // Get context ctx here (NVIDIA GeForce GTX 1080 Ti GPU) // ... // Copy data to GPU viennacl::compressed_matrix<double> A0(ctx); viennacl::copy(AA, A0); viennacl::vector<double> b0(b.size(), ctx); viennacl::copy(b, b0); // Instantiate solver viennacl::linalg::cg_tag my_cg_tag0(1e-5, 2000000); viennacl::linalg::cg_solver<viennacl::vector<double>> my_cg_solver0(my_cg_tag0); // Setup preconditioner viennacl::linalg::chow_patel_tag chow_patel_ilu_config; chow_patel_ilu_config.sweeps(3); chow_patel_ilu_config.jacobi_iters(2); viennacl::linalg::chow_patel_icc_precond<viennacl::compressed_matrix<double>> preconditioner(A0, chow_patel_ilu_config); // Solve system auto v0 = my_cg_solver0(A0, b0, preconditioner); Many thanks in advance, Rick |
From: Karl R. <ru...@iu...> - 2020-02-05 04:41:11
|
Hi Rick, have you verified that your matrices are positive definite? The problems with 'nan' usually stem from the lack of positive definiteness or zeros on the diagonal (as they often show up in saddle point problems). The code snippet you provide looks fine. 2000000 are a lot of iterations, though ;-) Best regards, Karli On 2/4/20 9:53 PM, Rick Fenrich wrote: > Hello all, > > I am solving large sparse symmetric systems. As a result I have been > using the CG solver. However the only preconditioner that seems to work > robustly across different matrices is the row_scaling preconditioner. I > have tried jacobi_precond and chow_patel_icc_precond with the CG solver > and while they work for some symmetric problems, they either fail to > converge or give NaNs for other symmetric problems. The same problem > converges just fine without the preconditioners or when using the > row_scaling preconditioner. > > Are there additional special settings beyond what is provided in the > documentation needed to successfully use these preconditioners? Are > there other limitations of these preconditioners? > > My matrices regularly contain non-zero entries O(10^6). I have verified > they are symmetric. I am using an NVIDIA GeForce GTX 1080 Ti GPU and > ViennaCL 1.7.1. > > Exemplary code is: > vector<map<unsigned int, double>> AA; /* Sparse matrix. */ > vector<double> b; /* Right hand side vector. */ > // Assemble AA and b here > // ... > // Get context ctx here (NVIDIA GeForce GTX 1080 Ti GPU) > // ... > // Copy data to GPU > viennacl::compressed_matrix<double> A0(ctx); > viennacl::copy(AA, A0); > viennacl::vector<double> b0(b.size(), ctx); > viennacl::copy(b, b0); > // Instantiate solver > viennacl::linalg::cg_tag my_cg_tag0(1e-5, 2000000); > viennacl::linalg::cg_solver<viennacl::vector<double>> > my_cg_solver0(my_cg_tag0); > // Setup preconditioner > viennacl::linalg::chow_patel_tag chow_patel_ilu_config; > chow_patel_ilu_config.sweeps(3); > chow_patel_ilu_config.jacobi_iters(2); > viennacl::linalg::chow_patel_icc_precond<viennacl::compressed_matrix<double>> > preconditioner(A0, chow_patel_ilu_config); > // Solve system > auto v0 = my_cg_solver0(A0, b0, preconditioner); > > Many thanks in advance, > > Rick > > > _______________________________________________ > ViennaCL-support mailing list > Vie...@li... > https://lists.sourceforge.net/lists/listinfo/viennacl-support > |
From: Rick F. <ri...@ar...> - 2020-02-05 22:44:21
|
Thank you for the quick reply Karl, The symmetric matrices are almost certainly positive definite because the system can be solved using the CG solver without preconditioning as well as using the CG solver with the row_scaling preconditioner. It is only the jacobi and chow_patel_icc preconditioners that cause difficulties. I have verified that there are no zeros on the matrix diagonal, but have not performed a positive definite check other than what I stated above. Thanks again for looking over the code snippet! All the best, Rick |
From: Rick F. <ri...@ar...> - 2020-02-07 23:03:42
|
After looking at the original paper on which the chow_patel_icc preconditioner is based, it seems that convergence for the preconditioner is guaranteed theoretically if a fixed point exists for the matrix. However, overflow can occur in practice which will lead to lack of convergence. In lieu of more detailed investigation, it seems likely that my matrix is inducing overflow during the preconditioner calculation which leads to lack NaNs during the CG solve. If any experts could comment on this possibility, please do so. Have others encountered this issue? |