Thanks again Matthew. I'll try to use MatSetNullSpaceFunction and see if it resolves the issue.
On Thu, Jun 1, 2017 at 6:00 PM, Matthew Knepley <[email protected]> wrote: > On Thu, Jun 1, 2017 at 4:58 PM, Bikash Kanungo <[email protected]> wrote: > >> Thank you Matthew for the quick response. >> >> I should provide some more details. >> >> 1) This implies that A is symmetric >> >> >> Yes A is symmetric >> >> 2) I think you mean the b is orthogonal to Q >> >> >> No. I've an extra condition on the solution x - that it has to be >> orthogonal to Q. Without this condition the problem will have multiple >> solutions. It's same as finding the minimum norm solution >> >> >>> Its not MatSetNullSpace() that fixes this issue, it is the linear solver >>> itself. Many Krylov solvers can converge to the >> >> minimum norm solution of this rank deficient problem. >> >> >> I was unable to get convergence without using MatSetNullSpace(). So it >> seems that the Krlov solvers are inadequate just by themselves. >> >> Second, we remove components in the nullspace from each iterate, which >>> is not what you are doing above. >>> >> >> This might be the reason why without MatSetNullSpace() I'm unable to >> attain convergence. >> >> >>> It seems easier to just give your Q as the nullspace. >> >> Yeah using Q as nullspace is the easiest option. But my basis is >> non-orthogonal and has an overlap matrix S. So in order to provide >> orthonormal nullvectors as Q to Petsc, I need to evaluate S^{-1/2}. As of >> now, for small problems I can evalaute S^{-1/2}, however, it's a show >> stopper for large ones. I can possibly use MatNullSpaceSetFunction() to >> circumvent the requirement of orthonormal nullvectors in MatSetNullSpace(). >> I would appreciate your input on this. >> > > You do not need this. Just use > > http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/Mat/ > MatNullSpaceSetFunction.html#MatNullSpaceSetFunction > > to implement your (I - Q^T Q). > > Thanks, > > Matt > > >> Lastly, I would like to know what operations do MatSetNullSpace enforce >> within a KSP solve. >> >> Thanks, >> Bikash >> >> On Thu, Jun 1, 2017 at 5:25 PM, Matthew Knepley <[email protected]> >> wrote: >> >>> On Thu, Jun 1, 2017 at 2:07 PM, Bikash Kanungo <[email protected]> wrote: >>> >>>> 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) This implies that A is symmetric >>> >>> 2) I think you mean the b is orthogonal to Q >>> >>> >>>> >>>> 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. >>>> >>> >>> Its not MatSetNullSpace() that fixes this issue, it is the linear solver >>> itself. Many Krylov solvers can converge to the >>> minimum norm solution of this rank deficient problem. Second, we remove >>> components in the nullspace from each >>> iterate, which is not what you are doing above. It seems easier to just >>> give your Q as the nullspace. >>> >>> Thanks, >>> >>> Matt >>> >>> >>>> Thanks, >>>> Bikash >>>> >>>> -- >>>> Bikash S. Kanungo >>>> PhD Student >>>> Computational Materials Physics Group >>>> Mechanical Engineering >>>> University of Michigan >>>> >>>> >>> >>> >>> -- >>> What most experimenters take for granted before they begin their >>> experiments is infinitely more interesting than any results to which their >>> experiments lead. >>> -- Norbert Wiener >>> >>> http://www.caam.rice.edu/~mk51/ >>> >> >> >> >> -- >> Bikash S. Kanungo >> PhD Student >> Computational Materials Physics Group >> Mechanical Engineering >> University of Michigan >> >> > > > -- > What most experimenters take for granted before they begin their > experiments is infinitely more interesting than any results to which their > experiments lead. > -- Norbert Wiener > > http://www.caam.rice.edu/~mk51/ > -- Bikash S. Kanungo PhD Student Computational Materials Physics Group Mechanical Engineering University of Michigan
