Boost logo

Ublas :

Subject: Re: [ublas] todo list
From: Matwey V. Kornilov (matwey.kornilov_at_[hidden])
Date: 2010-07-24 15:57:39

Nasos Iliopoulos wrote:

> To recapitulate my view is that we should provide custom implementations,
> open the development as possible so that people can contribute (maybe make
> an easy to access branch for all), and most importantly first define the
> interface it should have. The interface choices that come into my mind are
> the following (I will take SVD as an example): 1. Procedural:
> svd(A,U,S,V); 2. Functional - procedural: S = svd(A,U,V); 3. tieable:
> tie(U, S, V) = svd(A); tie(S) = svd(A); //or just S=svd(A); 4.
> Classical/or stl like: SVD(A.begin<1>(), A.end<1>(), S.begin<1>(),
> U.begin<1>(), V.begin<1>()) (ugly uh? fairly generic though).
> Overall I would vote for (3), but we need first to implement the "tie"
> interface.
I can't understand third option because a function can't be overloaded by
returned type.

However, there is one more well-known option:

void svd( Q, A, W ); // A := S, Q := U^T Q, W := W V

This is for memory economy. When we want to solve linear set Ax=y by means
of SVD decomposition we write

A=USV^T =>
USV^Tx=y =>
SV^Tx=U^Ty =>
x=V S^{-1} U^T y

and we see that (U^T y) can be calculated from y at place instead of
calculation U and performing matrix multiplication. Imagine that y is big
enough (say, 100), consequently U should be 100x100.