Hi, I'm trying to solve a linear system of equations Ax=b, where A has a null space (say Q) and x is known to be orthogonal to Q. In order to avoid ill-conditioning, I was trying to do the following:
1. Create A as a shell matrix 2. Overload the MATOP_MULT operation for with my own function which returns y = A*(I - QQ^T)x instead of y = Ax 3. Upon convergence, solution = (I-QQ^T)x instead of x. However, I realized that the linear solver can make x have any arbitrary component along Q and still y = A*(I-QQ^T)x will remain unaffected, and hence can cause convergence issues. Indeed, I saw such convergence problems. What fixed the problem was using MatSetNullSpace for A with Q as the nullspace, in addition to the above three steps. So my question is what exactly is MatSetNullSpace doing? And since the full A information is not present and A is only accessed through MAT_OP_MULT, I'm confused as how MatSetNullSpace might be fixing the convergence issue. Thanks, Bikash -- Bikash S. Kanungo PhD Student Computational Materials Physics Group Mechanical Engineering University of Michigan
