svd_lapack

From: Takeshi Enomoto (eno AT jamstec.go.jp)
Date: Mon Oct 24 2005 - 05:43:09 MDT


Hello,

I have a question on svd_lapack().

SVD is defined as,
      A = US(V^T)
so
      AV = US.

According to the NCL documentation,
V^T is returned when optv = 0 and V when optv = 1.

To confirm this I tested this function with the following simple example
to find two problems.

1. AV /= VS
2. A /= US(V^T) (but with optv = 1, A^T is obtained.)

Would you kindly point out where I made mistakes?

Takeshi

-------

begin
     n = 2
     a1d = (/.96, 1.72, 2.28, .96/)
     a = onedtond(a1d,(/n,n/))
     u = new((/n,n/), typeof(a))
     v = new((/n,n/), typeof(a))

     abuf = a
     optv = 1
     s = svd_lapack(abuf, "S", "S", optv, u, v)

     opt = True
     opt@title = "a"
     opt@row = True
     write_matrix(a, "2f7.2", opt)
     opt@title = "u"
     write_matrix(u, "2f7.2", opt)
     opt@title = "v"
     write_matrix(v, "2f7.2", opt)
     av = a#v
     opt@title = "AV"
     write_matrix(av, "2f7.2", opt)
     ss = (/(/s(0),0./),(/0.,s(1)/)/)
     opt@title = "S"
     write_matrix(ss, "2f7.2", opt)
     us = v#ss
     opt@title = "US"
     write_matrix(us, "2f7.2", opt)

     optv = 0
     vt = v
     abuf = a
     s = svd_lapack(abuf, "S", "S", optv, u, vt)
     aa = u#ss#vt
     opt@title = "A"
     write_matrix(aa, "2f7.2", opt)
end

_______________________________________________
ncl-talk mailing list
ncl-talk@ucar.edu
http://mailman.ucar.edu/mailman/listinfo/ncl-talk



This archive was generated by hypermail 2b29 : Wed Oct 26 2005 - 06:49:08 MDT