| Author |
Solving equations with multiple RHS
|
|
| Deepak 2008-02-28, 8:40 pm |
| Suppose I want to solve Ax = b_i, where A is a 500000x500000
sparse matrix for 1000 different RHS vectors (i.e., i
ranges from 1 to 1000). Something like x = A\b works. But is
this the best way to do this operation? How best can we do
this?
Thanks!
| |
| Bruno Luong 2008-02-28, 8:40 pm |
| "Deepak " <deepak.trivedi@gmail.com> wrote in message
<fq72dk$3jc$1@fred.mathworks.com>...
> But is
> this the best way to do this operation? How best can we do
> this?
I have no idea what "best" means. But A\B it's certainly a
very good way to treat multiple RHS.
Bruno
| |
| Deepak 2008-02-28, 8:40 pm |
| I tested it, and in fact it's amazingly fast. I'd be
interested in knowing what algorithm is used by Matlab for
doing multiple RHS.
Thanks.
"Bruno Luong" <b.luong@fogale.fr> wrote in message
<fq7597$j5i$1@fred.mathworks.com>...
> "Deepak " <deepak.trivedi@gmail.com> wrote in message
> <fq72dk$3jc$1@fred.mathworks.com>...
>
> I have no idea what "best" means. But A\B it's certainly a
> very good way to treat multiple RHS.
>
> Bruno
| |
| Petr Krysl 2008-02-28, 8:40 pm |
| "Bruno Luong" <b.luong@fogale.fr> wrote in message
<fq7597$j5i$1@fred.mathworks.com>...
> "Deepak " <deepak.trivedi@gmail.com> wrote in message
> <fq72dk$3jc$1@fred.mathworks.com>...
>
> I have no idea what "best" means. But A\B it's certainly a
> very good way to treat multiple RHS.
>
> Bruno
I'm not sure I agree. Yes, if A was attained by factoring
the original matrix; not necessarily otherwise.
Refer to a recent exchange on this topic:
http://www.mathworks.com/matlabcent...d/163743#417620
Petr
| |
| Tim Davis 2008-02-28, 8:40 pm |
| "Deepak " <deepak.trivedi@gmail.com> wrote in message
<fq7jdj$7je$1@fred.mathworks.com>...
> I tested it, and in fact it's amazingly fast. I'd be
> interested in knowing what algorithm is used by Matlab for
> doing multiple RHS.
>
> Thanks.
The matrix gets factorized once, then what happens next
depends on A and B, but the factors are used to do the
forward/backsolve.
Try:
spparms ('spumoni',2)
x=A\B ;
and it will tell you what it's doing.
| |
| Bruno Luong 2008-02-29, 4:41 am |
| "Petr Krysl" <pkryslNOSP@Mucsd.edu> wrote in message
<fq7ksj$mvs$1@fred.mathworks.com>...
>
> I'm not sure I agree. Yes, if A was attained by factoring
> the original matrix; not necessarily otherwise.
>
> Refer to a recent exchange on this topic:
>
>
http://www.mathworks.com/matlabcent...d/163743#417620
>
Yes, I have read the thread, not in deep though. It seems to
me the problem arises when transposing a triangular matrix
for back-solving, and this transposing requires memory and
time. OP asks here i the case where he/she has the A matrix
already built.
Or do I miss something else in the thread?
Cheers,
Bruno
| |
| Steven Lord 2008-02-29, 7:34 pm |
|
"Bruno Luong" <b.luong@fogale.fr> wrote in message
news:fq8gd3$sqk$1@fred.mathworks.com...
> "Petr Krysl" <pkryslNOSP@Mucsd.edu> wrote in message
> <fq7ksj$mvs$1@fred.mathworks.com>...
>
> http://www.mathworks.com/matlabcent...d/163743#417620
>
> Yes, I have read the thread, not in deep though. It seems to
> me the problem arises when transposing a triangular matrix
> for back-solving,
Why transpose? See step 3 in the Algorithms section of the MLDIVIDE and
MRDIVIDE reference page:
http://www.mathworks.com/access/hel...f/mrdivide.html
--
Steve Lord
slord@mathworks.com
| |
| Petr Krysl 2008-02-29, 7:34 pm |
| "Steven Lord" <slord@mathworks.com> wrote in message
<fq97bj$mhj$1@fred.mathworks.com>...
>
> "Bruno Luong" <b.luong@fogale.fr> wrote in message
> news:fq8gd3$sqk$1@fred.mathworks.com...
http://www.mathworks.com/matlabcent...d/163743#417620[color=darkred]
>
> Why transpose? See step 3 in the Algorithms section of
the MLDIVIDE and
> MRDIVIDE reference page:
Steve,
Now the problem there was that chol () would give you the
lower triangular factor, which then it needs to be used in
both the forward and backward substitution. However, from
my experiments it appears that to use the backward
substitution with the lower triangular matrix is not
possible without Matlab performing the transpose.
The key for an efficient solve would be to use
c=A\V; x = (c'/A)';
instead of
x = A'\(A\V);
Unfortunately, (c'/A) appears to be implemented as A'\c!
Since you volunteered your opinion here, which also be
willing to find out if this is the truth by looking under
the hood of the Matlab mldivide/mrdivide?
Thanks,
Petr
|
|
|
|