Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices |

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]

*To*: eigen@xxxxxxxxxxxxxxxxxxx*Subject*: Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices*From*: Benoit Jacob <jacob.benoit.1@xxxxxxxxx>*Date*: Wed, 24 Feb 2010 09:40:20 -0500*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=rgsxvGvczVTmbIHJF3y17P7iyRWD8i/2Q4RODcZISgI=; b=mSzxSUoPsycgwSBDualxgFPq6jEHbNLM3icY8oWYGn9wovGfxp/3/kVcecGEjkwhZm pv7z9tlX9ijVAxyobJYCxW+cYetopGHVNg9jgRQtI7HDIc6hX6F1M15PYdz9VvHw+Veg xeSvCRUj9e06dYS3d0cSp8O/l5YGcJ6vU2e98=*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=ptmqGhS7FzRuUTAXN/Tq0O+cfXtE0W1u+xES3wVJgD6Cpc89IJJHxG+Nh3IIw4VkQm KDoihj0MYPChIJRFZJWeyXhIpAdgJxS+quisO43uizjqoFNWspF5OUq/WED/UfSYKaG0 p4Or6nkB8gdhbeKYtc+PAy4+zzgWaq+kGSlMo=

2010/2/24 FMDSPAM <fmdspam@xxxxxxxxx>: > Thank you Benoit, > > (Sorry for my wired sentences in the previous mail. please forget the noise, > here is my second attempt) > > Am 24.02.2010 14:26, schrieb Benoit Jacob: >> >> 2010/2/24 FMDSPAM<fmdspam@xxxxxxxxx>: >> >>> >>> Hi Experts, >>> >>> please, for my personal clarification, let me ask: >>> (I'm not asking for a brief lesson in lin alg. but I have to happily live >>> with those answers :-) ) >>> >>> Def.: >>> D_labil := "how many zeros in vectorD()" >>> D_instable : "how many entries in vectorD() are smaller or equal 0.0" >>> >>> a) LDLt for non singular matrix is rank-revealing (of course, i hope), >>> >> >> This is a tautology: if you assume your matrices to be nonsingular, >> then you know their rank. Nonsingular means rank N. >> >> > > Yes, of course. It was listed for the sake of completeness. >>> >>> n_labil = 0 >>> >> >> Yes, for (sufficiently) nonsingular matrices, there won't be too small >> entries in the vectorD. >> >> >>> >>> b) LDLt for positive definite matrix => D_instable = 0 >>> >> >> Yes. >> >> >>> >>> c) LDLt for singular matrix is not rank-revealing, but: n_deficit> 0 >>> >> >> What is n_deficit ? Number of zeros in vectorD? > > Yes exactly, sorry. >> >> I think this kind of >> works, yes, as you're only interested in knowing if it's 0 or>0. As >> the LDLt decomposition still gives you the determinant, it in >> particular tells you very vaguely if the matrix is singular or not, >> but don't use it for that, it shouldn't be very precise. >> >> > > A better alternative exist, right? (see my (ab-)use case below). If you want to measure the rank (deficit), you need a rank-revealing decomposition. For example: - full-pivoting LU - column-pivoting QR (yes it's enough) - full-pivoting QR - SVD >>> >>> d) LDLt for indefinite matrix => D_instable> 0 but counting D(i)<= 0.0 >>> is >>> meaning less? >>> >> >> Yes. The "zeros" in the vectorD will actually be small positive >> values. There's no magic cutoff where you'd decide that a small value >> can be considered zero. >> >> > > OK. Maybe from a practical point of view: D(i,i) < D.max*1e-8 could > considered as "small value"? But that may depend on concrete use cases. I don't know, this could be trickier depending on how big the off-diagonal elements of the L matrix can get. All I can do is point you to the above Higham article, I'm not an expert myself. > My (ab-)use case or "how semi-education cause numerical nightmares": > > I have infiltrated a very ancient but nevertheless successfully just in use > Fortran FE-code with my own taucs-LLt-solver binding (with a very welcome > performance boost). > The taucs-LLt-solver itself use blockwise a traditional Lapack-solving > routine. > > In that Lapack solving routine, I do observe and count "nearly" negative or > fully diagonal elements (A(i,i) < +eps). > If I detect a "negative" diagonal elements, I increment my a counter (say > "++D_instabile") and brutally set the diagonal element to A(i,i) = +1.0.. the > rest of the Lapack routine thereon is used without further modifications. > That count of "negative" diagonal elements is used as an indicator for an > outer iteration loop to find the instability limit (D_instabile = 0 <-> > D_instabile=1). > > This algorithm is in use, but tend to be unstable for bigger number of > degrees of freedom. > Now, thank you, I have several new staring points. Good luck :) Benoit > > Frank > > >

**References**:**[eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Ben Goodrich

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Ben Goodrich

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Benoit Jacob

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Ben Goodrich

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Gael Guennebaud

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Jitse Niesen

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Benoit Jacob

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*FMDSPAM

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*Benoit Jacob

**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices***From:*FMDSPAM

**Messages sorted by:**[ date | thread ]- Prev by Date:
**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices** - Next by Date:
**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices** - Previous by thread:
**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices** - Next by thread:
**Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices**

Mail converted by MHonArc 2.6.19+ | http://listengine.tuxfamily.org/ |