|Re: [eigen] Instability in LLT and LDLT methods.|
[ Thread Index |
| More lists.tuxfamily.org/eigen Archives
- To: eigen@xxxxxxxxxxxxxxxxxxx
- Subject: Re: [eigen] Instability in LLT and LDLT methods.
- From: Benoit Jacob <jacob.benoit.1@xxxxxxxxx>
- Date: Wed, 28 Jan 2009 16:34:16 +0100
- Dkim-signature: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=domainkey-signature:mime-version:received:in-reply-to:references :date:message-id:subject:from:to:content-type :content-transfer-encoding; bh=LbN+bTXC7KQPlZABoeWM7l+0zp9vt87mLdjlqsC5q8o=; b=ZnKBmmKhGlzdqkdux6NCazUy/rFChycHXQEE2+ma/PLqFtxvjG8N8jms7gD7J2yOTv QNyfG0794B3MpwBX+3OKEcqepvlKK09c0JBppJ5hTt3vUeWEOzCVI+w6IVniv10cr+pZ Z8josj0LpNgdz9Wv0mYMg5raBOsN2Nx0lmphs=
- Domainkey-signature: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=mime-version:in-reply-to:references:date:message-id:subject:from:to :content-type:content-transfer-encoding; b=Dw9Lp3MJf+Fmq5PTGF8iBDzR/Iu9VR1ZCZS5pgy0Pfss8ZDVU+j5Ezp1BO5GEzzD3M QWLYuDCCQqN/f3YZ3Q3bB0z/4mdGCATsSIRXRUbuFREGneEvdKEv1ZfxEtcuNkxyK3TV tOE1BnyaIhoteS+b2EecA2Flgc875CDZXUqjk=
2009/1/28 Jitse Niesen <jitse@xxxxxxxxxxxxxxxxx>:
> On Wed, 28 Jan 2009, Benoit Jacob wrote:
>>> When solving a linear system with an invertible matrix, partial pivoting
>>> performs as well as complete pivoting, doesn't it? LAPACK uses partial
>>> pivoting, as does GSL.
>> No, partial pivoting fails even for invertible matrices of size
>> 200x200, and fails consistently at sizes like 1000x1000.
> How did you arrive at this conclusion?
When I implemented LU, I was initially hoping that partial pivoting
would be enough. What I say here is the results of my experiments. So
maybe I got something wrong, but aside from that possibility, I was
talking first-hand experience when saying that.
>I'm by no means an expert in
> numerical linear algebra, but the books I have (including Golub & Van Loan,
> Matrix Computations, 3e and Higham, Accuracy and Stability ..., 2e) say that
> partial pivoting should work fine.
At least Golub & van Loan do mention the numerical stability problem
of partial pivoting...
(plus I don't consider this book to be a solid reference on numerical
stability, e.g. with Keir we found a numerical stability issue in
their householder vector computation).
>> I think that LAPACK and GSL do a hybrid of partial and full pivoting
>> that I have yet to understand.
> That's not what their documentation says:
Ah OK. Here I was just repeating stuff I heared over IRC. So I don't
know. Perhaps they still do kind of a hybrid but don't mention it in
the docs. If they really do plain partial pivoting, then let's try
inverting, say, a 2000x2000 matrix, it will fail. For this reason I
doubt that they actually do plain partial pivoting.
> I think a 50% speedup is a worthwhile target. At least for vector solves,
> you can compute the residual and redo the computation with complete pivoting
> if necessary.
Here we're talking about a cubic complexity algorithm. So while a 50%
speedup sounds big, it just means that in the same time you handle a
matrix with 14% more rows and columns ( 1.5^(1/3) = 1.14). I want
Eigen to be easy to use and that implies prefering safeness over
Now let me retry my experiment with partial pivoting to see if it's