[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

*Subject*: Re: Cholesky factorisation*From*: "David McClain" <dmcclain(at)azstarnet.com>*Date*: Tue, 22 Feb 2000 00:28:59 -0700*Newsgroups*: comp.lang.idl-pvwave*Organization*: Posted via Supernews, http://www.supernews.com*References*: <01bf7cda$4593eba0$eb2f14ac@default>*Xref*: news.doit.wisc.edu comp.lang.idl-pvwave:18529

I use the method shown in Numerical Recipes, 2nd Ed. Coding in C and being clever allows it to run in under 3 microsec on a 6x6 array, which is a size commonly needed by one of our analyses. Doing it externally allows for both speed and error handling as you would want it. It also guarantees correctness, as verified by you, the coder; something that RSI would have you accept as a matter of faith... David McClain, Sr. Scientist Raytheon Systems Co. Tucson, AZ Peter Scarth <p.scarth@mailbox.uq.edu.au> wrote in message 01bf7cda$4593eba0$eb2f14ac@default">news:01bf7cda$4593eba0$eb2f14ac@default... > Howdy all, > I'm trying to determine if a symmetric matrix is positive-definite. If this > sounds like gobbdlygook you might like to stop reading here. If the > idl2matlab translate-o-matic existed, it would probably need to know about > this to translate the handy \ (mldivide) operator. There a number of tests > for positive-definiteness and the one I'm looking at is to attempt cholesky > factorisation of the matrix. IDL's CHOLDC procedure halts and returns > CHOLDC: choldc failed (note that matlab can return a value to flag a > failure in the decomposition). > The metho that I an using to get the decomposition looks something like > this: > > -------------- > pro choldc2,a,p > sa=size(a) > nd=sa[0] & m=sa[1] & n=sa[2] > > ; Some initial tests..... > ; ..... > > p=fltarr(n) > for i=0,n-1 do begin > for j=i,n-1 do begin > sum=a(j,i) > if i gt 0 then sum=sum-total(a(0:i-1,i)*a(0:i-1,j)) > if (i eq j) then begin > if (sum le 0) then begin > p=-1 > return > endif else p(i)=sqrt(sum) > endif else a(i,j)=sum/p(i) > endfor > endfor > return > end > ---------------------- > > it works ok, returns -1 in p if the method fails but only runs at 1/10 the > speed of the built in version. > Is it possible to vectorise this further, or has someone out there in > cyberland found a more elegant solution to this problem already? > > Thanks, > > > Peter Scarth > Biophysical Remote Sensing Group > The University of Queensland. >

**Follow-Ups**:**Re: Cholesky factorisation***From:*Pavel Romashkin

**References**:**Cholesky factorisation***From:*Peter Scarth

- Prev by Date:
**Re: idl2matlab translate-o-matic** - Next by Date:
**Re: idl2matlab translate-o-matic** - Prev by thread:
**Cholesky factorisation** - Next by thread:
**Re: Cholesky factorisation** - Index(es):