One more question.
Suppose that we are solving the singular linear system Ax = b. N(A) is the
null space of A, and N(A^T) is the null space of the transpose of A.

## Advertising

The linear system is solved using SNES, that is, F(x) = Ax-b = Ax -b_r -
b_n. Here b_n in N(A^T), and b_r in R(A). During each nonlinear
iteration, a linear system A \delta x = F(x) is solved. N(A) is applied to
Krylov space during the linear iterating. Before the actual solve
"(*ksp->ops->solve)(ksp)" for \delta x, a temporary copy of F(x) is made,
F_tmp. N(A^T) is applied to F_tmp. We will get a \delta x. F(x+\delta x )
= A(x+\delta x)-b_r - b_n.
F(x+\delta x ) always contain the vector b_n, and then the algorithm never
converges because the normal of F is at least 1.
Should we apply N(A^T) to F instead of F_tmp so that b_n can be removed
from F?
MatGetTransposeNullSpace(pmat,&nullsp);
if (nullsp) {
VecDuplicate(ksp->vec_rhs,&btmp);
VecCopy(ksp->vec_rhs,btmp);
MatNullSpaceRemove(nullsp,btmp);
vec_rhs = ksp->vec_rhs;
ksp->vec_rhs = btmp;
}
should be changed to
MatGetTransposeNullSpace(pmat,&nullsp);
if (nullsp) {
MatNullSpaceRemove(nullsp,ksp->vec_rhs);
}
???
Or other solutions to this issue?
Fande Kong,
On Thu, Oct 13, 2016 at 8:23 AM, Matthew Knepley <knep...@gmail.com> wrote:
> On Thu, Oct 13, 2016 at 9:06 AM, Kong, Fande <fande.k...@inl.gov> wrote:
>
>>
>>
>> On Wed, Oct 12, 2016 at 10:21 PM, Jed Brown <j...@jedbrown.org> wrote:
>>
>>> Barry Smith <bsm...@mcs.anl.gov> writes:
>>> > I would make that a separate routine that the users would call first.
>>>
>>> We have VecMDot and VecMAXPY. I would propose adding
>>>
>>> VecQR(PetscInt nvecs,Vec *vecs,PetscScalar *R);
>>>
>>> (where R can be NULL).
>>>
>>
>> What does R mean here?
>>
>
> It means the coefficients of the old basis vectors in the new basis.
>
> Matt
>
>
>> If nobody working on this, I will be going to take a try.
>>
>> Fande,
>>
>>
>>>
>>> Does anyone use the "Vecs" type?
>>>
>>
>>
>
>
> --
> 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
>