Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Solving large sparse non-square system

Subject: Solving large sparse non-square system

From: Nargol

Date: 24 Nov, 2010 21:07:05

Message: 1 of 5

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

Subject: Solving large sparse non-square system

From: Matt J

Date: 24 Nov, 2010 21:47:03

Message: 2 of 5

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
 

Subject: Solving large sparse non-square system

From: Bruno Luong

Date: 25 Nov, 2010 08:01:08

Message: 3 of 5

"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

Subject: Solving large sparse non-square system

From: Nargol

Date: 26 Nov, 2010 03:49:04

Message: 4 of 5

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

Cheers

Subject: Solving large sparse non-square system

From: Brian Borchers

Date: 26 Nov, 2010 05:06:13

Message: 5 of 5

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.

Tags for this Thread

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us