## Glas :## [glas] (no subject) |

**From:** Wolfgang Bangerth (*bangerth_at_[hidden]*)

**Date:** 2005-10-13 19:01:50

**Next message:**Laurent Plagne: "Re: [glas] (no subject)"**Previous message:**Neal Becker: "Re: [glas] norms etc"**Next in thread:**Laurent Plagne: "Re: [glas] (no subject)"**Reply:**Laurent Plagne: "Re: [glas] (no subject)"**Reply:**Karl Meerbergen: "Re: [glas] (no subject)"

In reply to an already several days old message:

> I recall a demand for supporting vectors of matrices and matrices of

> matrices. My opinion is that this should be another type of collection

> where the value_type is still the value_type of the lowest level of

> matrix in the nesting of matrices.

There is definitely a demand for something like this, though I wouldn't

suggest that implementing them as matrices of matrices is a particularly

clever way to go.

The application where this is important is when one can split a matrix

into blocks. Typical applications are discretizations of vector-valued

partial differential equations (think: mixed methods, saddle point

problems) where instead of solving with a matrix directly, one would like

to form Schur complements of its blocks. The right abstraction therefore

is to consider the whole matrix as an array of submatrices (blocks).

The way we implement this in deal.II (and I believe this is the correct

way) is indeed as an array of matrices. One can access individual blocks,

but if you try to access individual elements of the matrix, one gets a

reference to the corresponding element in the block in which it resides.

The BlockMatrix class therefore has a dual interface: access the blocks of

the array through

Matrix & block (const unsigned int block_i, const unsigned int block_j)

and access the elements of the matrix through something like

reference operator() (const unsigned int i, const unsigned int j) {

return blocks[block_from_global(i)][block_from_global(j)]

->operator() (local_from_global(i), local_from_global(j));

}

Through the second interface, one can hide the actual storage scheme from

algorithms that only expect a matrix, while the first interface allows to

extract and operate on individual blocks.

If, on the contrary, this would be implemented as matrices-of-matrices,

then one would only have the first interface, but algorithms that just

expect a transparent matrix object won't work because your operator()

doesn't return individual entries of the matrix, but an entire block.

For an example of an implementation of such concepts, feel free to take a

look at

http://www.dealii.org/5.2.0/doxygen/lac/classBlockSparseMatrix.html

Best

Wolfgang

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

Wolfgang Bangerth email: bangerth_at_[hidden]

www: http://www.math.tamu.edu/~bangerth/

**Next message:**Laurent Plagne: "Re: [glas] (no subject)"**Previous message:**Neal Becker: "Re: [glas] norms etc"**Next in thread:**Laurent Plagne: "Re: [glas] (no subject)"**Reply:**Laurent Plagne: "Re: [glas] (no subject)"**Reply:**Karl Meerbergen: "Re: [glas] (no subject)"