Thread Subject:

Solving large sparse non-square system

Hello,

I have a large sparse non-square linear system, I need to solve. Back slash operator(mldivide) does not help due to memory problems. I have also used lsqr, but it is does not return a good enough solution (back slash works better, but it is not useful for large problems). I know about the iterative methods such as

    pcg - Preconditioned Conjugate Gradients Method.
    bicg - BiConjugate Gradients Method.
    bicgstab - BiConjugate Gradients Stabilized Method.
    cgs - Conjugate Gradients Squared Method.
    gmres - Generalized Minimum Residual Method.
    minres - Minimum Residual Method.
    qmr - Quasi-Minimal Residual Method.
    symmlq - Symmetric LQ Method.

But non of them are appropriate for my problem, since my coefficient matrix is non-square. Could you please help me find an appropriate solver for my problem?

Thanks,

Nargol

If your given system is A*x=b, you could try solving the equivalent system

A'*A*x=A'*b;

which is a square system. It may not work so well because it will generally have poorer conditioning. The following alternative algorithm also sometimes works well:


absA=abs(A);
C=absA*sum(absA,2);
x=initialValues;

for ii=1:numIter

 x=x-A.'*(A*x-b)./C;

end
 

"Matt J " <mattjacREMOVE@THISieee.spam> wrote in message <ick14n$5b1$1@fred.mathworks.com>...
> If your given system is A*x=b, you could try solving the equivalent system
>
> A'*A*x=A'*b;
>
> which is a square system. It may not work so well because it will generally have poorer conditioning.

I see may people write such statement and to my view it is also half true, for two reasons:

If the matrix A has badly conditioning so that the square of that has numerical problem, then user play already with the fire. The original problem is already badly conditioned and it needs to be taken care by other mean (regularization, adding constraints, etc...), but then the condition number will reduces.

With an iterative methods such as conjugate gradient - I believe it is not matter to multiply by A' and stretch the spectrum, since the Krylov space contains the power (A'*A)^k*x and not the inverse (A'*A)^-k*x.

For large system, solving normal equation is sometime the only option. And it works better than you might expect. Many sophisticated weather forecast codes take into account the data using normal equation.

Bruno

Thanks for your reply. I will try that. Hope it works.

Cheers

On Nov 24, 2:07 pm, "Nargol " <x...@x.com> wrote:
> Hello,
>
> I have a large sparse non-square linear system, I need to solve. Back slash operator(mldivide) does not help due to memory problems. I have also used lsqr, but it is does not return a good enough solution (back slash works better, but it is not useful for large problems).

It would be worth while to understand how LSQR is failing. If LSQR is
converging, but the solution is not acceptable because ||Ax-b|| is too
large, then you have a more fundamental problem, because LSQR is
solving the least squares problem and you're just not satisfied with a
least squares solution even though there isn't an exact solution to
your problem.

On the other hand, if LSQR is failing to converge (maximum iterations
exceeded), than this is an indication that the least squares problem
is ill-conditioned. In that case, some kind of regularization would
typically be the appropriate thing to try.