|
|
Ublas : |
Subject: Re: [ublas] symmetric matrices, row and column-major
From: Tiago Requeijo (tiago.requeijo.dev_at_[hidden])
Date: 2009-02-24 17:08:36
Looking at functional.hpp, I think the problem is indeed the
transposed_structure template. As far as I understand, the template simply
takes an upper matrix and regards it as the transpose of a lower matrix.
While this is (of course) correct, it doesn't take into account the
row/column major storage:
If we take the upper matrix U
1 2 3
- 4 5
- - 6
and regard it as the transpose of a lower matrix L
1
2 4
3 5 6
without a row/column major swap, then we're storing the row major version of
the lower matrix or, in other words, the column major version of the upper
matrix:
a) row_major for U: 1 2 3 4 5 6
b) row_major for L: 1 2 4 3 5 6
Since the transposed_structure template simply maps the elements to the
transpose and doesn't change the storage, we are in fact storing b) instead
of a)
Tiago
On Tue, Feb 24, 2009 at 2:11 PM, Gunter Winkler <guwi17_at_[hidden]> wrote:
> Tiago Requeijo schrieb:
> > I'm having a problem with row and column major in symmetric matrices.
> > I'm probably not understanding correctly how the storage takes place
> > in memory. I assumed upper+column_major storage was the same as
> > lower+row_major. Here's an example:
>
> The mapping from (i,j) to storage location is defined in
> basic_row_major::lower_element and basic_row_major::upper_element (see
> functional.hpp).
>
> symmetric_matrix::operator(i,j) uses the triangular_type::element(...)
> method to dispatch to the corresponding storage mapping from the layout
> type (row major/column major). Maybe this mapping does not work as
> expected.
>
> I recently introduced the transposed_structure template to reduce the
> number of implementations for the diffentent triangular types. Can you
> check if the problem is caused by a wrong implementation of
> transposed_structure::element()?
>
> > Upper part of the symmetric matrix:
> >
> > 0.34 0.9 1.17 0.689 1.42
> > ---- 1.36 0.9 0.808 0.146
> > ---- ---- 0.572 1.55 1.61
> > ---- ---- ---- 0.093 1.32
> > ---- ---- ---- ---- 0.828
> >
> >
> > ublas::symmetric_matrix<double, ublas::upper, ublas::column_major> S
> > double* a = &(S.data()[0]);
> > ublas::symmetric_matrix<double, ublas::lower, ublas::column_major> S2
> > double* a2 = &(S2.data()[0]);
> >
> > Contents of a, a2:
> > 0.34 0.9 1.17 0.689 1.42 1.36 0.9 0.808
> > 0.146 0.572 1.55 1.61 0.093 1.32 0.828
> >
> This seems to be wrong ...
>
> mfg
> Gunter
> _______________________________________________
> ublas mailing list
> ublas_at_[hidden]
> http://lists.boost.org/mailman/listinfo.cgi/ublas
>