## Glas :## Re: [glas] multi-dimensional arrays? |

**From:** Peter Gottschling (*pgottsch_at_[hidden]*)

**Date:** 2006-12-18 09:28:48

**Next message:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Previous message:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"**In reply to:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"**Next in thread:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Reply:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Reply:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"

Hi Laurent,

Concerning nested matrix and vector types, we need to ask ourselves

first if we need special algorithms to handle it or if a least in some

cases the operators on the values already use the right algorithm.

Second, is the use of operators efficient enough performance-wise? In

a matrix multiplication c(i, k)+= a(i, j) * b(j, k) is correct for any

type with consistent operators but not very efficient if a, b, and c

are matrices. Unless we have really sophisticated expression templates.

Next question, do we need special concepts? Say e.g.

RecursiveTransposable for a recursive function transpose

How would it look like?

concept RecursiveTransposable<typename Matrix>

{

where RecursiveTransposable<Matrix::value_type>;

// ...

};

would recursively descent and break if it reaches a scalar value that

has no value_type.

Conversely, something like:

auto // ???

concept RecursiveTransposable<typename Matrix>

{

Matrix transpose(Matrix );

}

Now we say this is defined for arithmetic scalars

template<typename T>

where Arithmetic<T>

concept_map RecursiveTransposable<T> {}

And now enable it from bottom up:

template<typename T>

where RecursiveTransposable<T>

concept_map RecursiveTransposable<matrix<T> > {}

Extending this to block-matrices,

(A B)

(C D)

where A, B, C, and D are matrices of different type. Say we have a type

for heterogeneous blocks

template <typename MA, typename MB, typename MC, typename MD>

hb_matrix { ... };

We can define:

template <typename MA, typename MB, typename MC, typename MD>

where RecursiveTransposable<MA>

&& RecursiveTransposable<MB>

&& RecursiveTransposable<MC>

&& RecursiveTransposable<MD>

concept_map RecursiveTransposable<hb_matrix<MA, MB, MC, MD> > {}

However, this are only some preliminary thoughts. We need to look,

which algorithms we want to implement recursively. And it looks like

that there can be a great difference between a syntactically correct

implementation and a high-performance implementation.

Cheers,

Peter

On 15.12.2006, at 05:49, plagne.laurent_at_[hidden] wrote:

> Hi Karl,

>

> IMHO, the matrix< matrix<..> > can be more or less difficult to

> implement

> depending on the constraints you accept. As an example, if you decide

> that the

> library must remain open to third party codes and handle different

> kind of

> matrix data storage : this is pretty easy to do with simple matrix type

> (matrix<double>) because the return type of A(i,j) will be a simple

> ref to

> double (or const ref). In the 2D case, matrix< matrix<double> > A, the

> accessor

> operator A(i,j) will be naturally a ref to matrix<double>. And you

> loose the

> compatibility with third party data elements...

> One can avoid this problem but it is not a trivial issue.

>

> Another point in nested matrix issue is that it can help (again IMHO)

> to

> validate the very impressive work the GLAS team has done with vector

> space

> concepts (HI Peter). For example, the 0 and 1 matrix counter-part are

> to be

> efficiently treated.

>

> Last point is maybe not a so general issue : In my case, I have to

> handle sparse

> matrix with a complicated but highly structured sparsity pattern. This

> pattern

> can be express with a reduced set of simple sparsity patterns : e.g A

> =Dense<Sparse<Diagonal..> > > + Sparse<Tridiagonal<Dense<..> > >

>

>

> Best

>

> Laurent

>

>

>> I think this is an important point. Toon and I have thought about this

>> for a long time. We will certainly talk about this in detail.

>> Providing

>> the container matrix< matrix<T> > is not hard. I recall that

>> algorithms

>> are application/interpretation dependent, and we had a lively

>> discussion

>> on the list on this topic. So, I agree that 'some' support is needed,

>> although it is not yet clear to me what 'some' here means.

>>

>> Best,

>

>

>

> _______________________________________________

> glas mailing list

> glas_at_[hidden]

> http://lists.boost.org/mailman/listinfo.cgi/glas

>

------------

Peter Gottschling, Ph.D.

Research Associate

Open Systems Laboratory

Indiana University

135 Lindley Hall

Bloomington, IN 47405

Tel.: +1-812-855-3608 Fax: +1-812-856-0853

http://www.osl.iu.edu/~pgottsch

**Next message:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Previous message:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"**In reply to:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"**Next in thread:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Reply:**Karl Meerbergen: "Re: [glas] multi-dimensional arrays?"**Reply:**plagne.laurent_at_[hidden]: "Re: [glas] multi-dimensional arrays?"