Cholesky factorisation

["Peter Scarth" <p.scarth(at)mailbox.uq.edu.au>]

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.