home *** CD-ROM | disk | FTP | other *** search
/ NetNews Usenet Archive 1992 #20 / NN_1992_20.iso / spool / comp / theory / 1904 < prev    next >
Encoding:
Internet Message Format  |  1992-09-09  |  2.2 KB
  1. Path: sparky!uunet!mcsun!uknet!cam-cl!cam-cl!nmm
  2. From: nmm@cl.cam.ac.uk (Nick Maclaren)
  3. Newsgroups: comp.theory
  4. Subject: Re: Matrix Square Root
  5. Message-ID: <1992Sep10.075331.9873@cl.cam.ac.uk>
  6. Date: 10 Sep 92 07:53:31 GMT
  7. References: <Kevin.D.Brunson.1-080992181444@davmac27.oshaughnessy.lab.nd.edu> <dak.716054795@kaa>
  8. Sender: news@cl.cam.ac.uk (The news facility)
  9. Reply-To: nmm@cl.cam.ac.uk (Nick Maclaren)
  10. Organization: U of Cambridge Computer Lab, UK
  11. Lines: 41
  12.  
  13. In article <dak.716054795@kaa>, dak@kaa.informatik.rwth-aachen.de (David
  14. Kastrup) writes:
  15. > Kevin.D.Brunson.1@nd.edu (Kevin D. Brunson) writes:
  16. > > ....
  17. > >        The problem I have seems like a simple one:  How do I calculate the
  18. > >square root of a symetric matrix?  Theoretically, if a matrix A is positive
  19. > >semidefinite then there exists a matrix B that is the square root of A such
  20. > >that A = B*B.  I have consulted Lancaster and Tismenetsky's "The Theory of
  21. > >Matrices" but a programmable algorithm is not obvious.  ....
  22.   
  23. > If you do an Eigenvector decomposition A=TDT^{-1}, where D is a diagonal
  24. > matrix containg the eigenvalues of the matrix, and T a matrix containing
  25. > the corresponding eigenvectors in its columns, TD^{0.5}T^{-1} will
  26. > have that property.
  27.  
  28. I hit this problem when writing the multivariate normal random number
  29. generator for the NAG library, and there are three solutions:
  30.  
  31.     1) Gaussian elimination on the original matrix plus a fudge constant
  32. on the diagonals (which I vaguely remember as being 5*maximum
  33. element*epsilon).
  34.  
  35.     2) Cholesky decomposition ditto (but the fudge constant is larger,
  36. possible 13*...).  This is the method that the NAG code uses.
  37.  
  38.     3) A = T'DT decomposition, as described above, and replace negative
  39. elements of D by zero.
  40.  
  41.     The speed stakes are (2) > (1) > (3) and the accuracy stakes are
  42. exactly the converse.  In all cases, the difference is a small constant
  43. factor, so there is nothing much in it.  For the theory of (1) and (2),
  44. see Wilkinson and Reinsch - I don't know much about (3).
  45.  
  46.  
  47. Nick Maclaren
  48. University of Cambridge Computer Laboratory,
  49. New Museums Site, Pembroke Street,
  50. Cambridge CB2 3QG, England.
  51. Email:  nmm@cl.cam.ac.uk
  52. Tel.:   +44 223 334761
  53. Fax:    +44 223 334679
  54.