Re: ublas: hermitean matrices

From: walter_at_[hidden]
Date: 2001-12-02 04:39:35

--- In boost_at_y..., Peter Schmitteckert (boost) <boost_at_s...> wrote:

[snip]

> > >   I'd prefer the diagonal elements of a complex hermitean
> > >   matrix to be real. Or at least to have a function
> > >   returning the diagonal elements as real values.
> >
> > How would you design the corresponding iterator interface using 
real
> > diagonal elements? Shouldn't we better assert, that imag (m [i] 
[i])
> > == 0? I'll think about an special accessor for diagonal elements.
>
> In my own matrix package I did the following:
>    1) a special function that return diagonal elements, e.g. return 
type 
> double for complex<double>
>   2) the element access through
>     F& HermitianMatrix::operator(const i, const j)
>     is only defined for i<j

I'm rather sure, that we should try to go with double d = real (hm 
(i, i)). Isn't that similar enough to double d = hm.diag (i)? There 
also is an open issue with proxies like diag, which I'll address in 
another posting.
 
>   but I had no iterators defined :(, so I think a special accessor
>   for diagonal elements is the way to go.
> 
> > > -> In the same spirit one should not be allowed to
> > >   multiply a complex hermitean matrix with a non-real
> > >   scalar, or at least this shouldn't return a hermitean matrix.
> >
> > The most critical part here is hermitean_matrix *= scalar. Again 
I'd
> > prefer to assert that imag (scalar) == 0.
>
> But that's sensitive to rounding problems, e.g. if you're writing 
an eigen
> system solver you might create minimal imaginary parts, and you have
> to set the explicitly to zero.

The type of the scalar in hermitean_matrix *= scalar is a template 
argument. So it should be no problem to write hm *= d;
 
> > I imagine, that deriving banded_hermitean_matrix (or
> > hermitean_banded_matrix :-) statically from matrix_expression 
should
> > give you the neccessary ET support.
> >
> > A hint: may be you better try to implement a somewhat simpler form
> > like skew_hermitean_matrix ;-) first.
> >
> > Thanks for your feedback.
> 
> o.k., maybe you're right, 
> would you be willing to audit/comment my attempts to implement
> new matrix types ?

Yes (although we're losing one of our rare users ;-).
 
Regards
 
Joerg

P.S.: Did you consider to use hermitean_adaptor<banded_matrix<> >?