Re: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type? |

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]

*To*: eigen <eigen@xxxxxxxxxxxxxxxxxxx>*Subject*: Re: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?*From*: Gael Guennebaud <gael.guennebaud@xxxxxxxxx>*Date*: Sun, 17 Mar 2013 19:38:13 +0100*Dkim-signature*: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=20120113; h=x-received:mime-version:in-reply-to:references:from:date:message-id :subject:to:content-type:content-transfer-encoding; bh=Un+/JFnMyi6aE9UpgA9jbf6IXDiOlHbcHxmkenz4wb4=; b=meSywvEvdFkhsoHQlhFP2nnxfdsW8xC/gMlpEEhQ4iDI3UD5doOXC2ujA5IG2lZlQs 0QOwGjIS4IBZk5MVkKdNEmBOOf0c/+OXjzVD4JdRnzaYmbUCQeLa948VPIN74Jfxl+8b 5wO1EwV/f42Sv673oKb/43YIEHfIVD7NHQncDrlxR0j2LUafZczvWtkM0PJ1zNFmLvYR 5KWwQwj/7vTxEmOa7yas36cfB5DaJw3dw8EQe/q6XdKnxdYHWuJC7WsjFbgLPVk5JC3U qPlKkIYlvKvzi+aJzZe1BKuhT7Ie+CqMRzwHG1OrYswEXm8e/ebHCe9Zo4AosNLfaZ6p rmYA==

Hi Daniel, On Sat, Mar 16, 2013 at 7:00 PM, <Daniel.Vollmer@xxxxxx> wrote: > I've attached my experiment that essentially uses a dense 5x5 matrix as the "Scalar" type. Basically, this works (which is already pretty nice) but I have noticed two things which seem sub-optimal, and so I was curious whether I'm missing something or doing it wrong. adding an explicit support for sparse block matrices is something we planed to address soon. Your example is nice starting proof-of-concept. To get it right, we certainly need a new SparseBlockMatrix class and some way to specify which algorithm to use for a/b, i.e., a dense Cholesky in SimplicialCholesky, LU in SparseLU, etc. > 1) The Scalar-type used in SparseMatrix seems to need an abs() function in order to (I assume) only iterate over non-zero elements (AmbiVector.h:308 & 335). Oddly enough, my abs-version was needed to compile, but never got called. > It might be advantageous to provide an iterator that iterates over all elements without making such a guarantee, but then the point shouldn't matter once it is compressed. abs is probably not the only free fonction you/we will have to specify, but also operator/, and in some cases sqrt (LLT). In this precise case of AmbiVector, is only needed for non-conservative sparse matrix products, i.e.: (A*B).prune(eps); > 2) Many times the single-argument constructor is called for the Scalar-type of SparseMatrix in instances where I don't think this is necessary. As this becomes more and more expensive with a more heavy-weight Scalar type, I'd like to avoid this where possible, for example > - SparseMatrix.h:1108 (insertUncompressed): Element is initialised to 0 even though it is later overwritten by the inserted value. Worked fine for me if I removed the = 0. > - SparseMatrix.h:382 (insertBackByOuterInner) and SparseMatrix.h:392 (insertBackByOuterInnerUnordered): Element is initialised to 0 even though it is later overwritten by the inserted value. I hacked this by making a variant of CompressedStorage::append() that only takes an Index and not a value-ref and calling that instead (because again, the value is later overwritten via the reference we return) > - AmbiVector.h:197/206/239 (coeffRef): I don't quite understand this, but I'm not sure the element needs to be initialised here either. We need these initializations because it is allowed to do: mat.insert(i,j); to insert an explicit zero. One solution could be to return a wrapper object that would set the wrapped value to zero in the destructor if operator= has not been called, but that sounds a bit overkill. > - AmbiVector.h:171 (setZero): Shouldn't / couldn't this use assignment from m_zero instead? Actually, all these constant-to-Scalar conversions are problematic for any heavy-weight scalar types (AutoDiff, multi-precision, etc.), and an old idea was to add in NumTraits a Constant (maybe that name is not good) type for that usage. This way if custom scalar type is fully compatible with double, then NumTraits<Custom>::Constant could be double. For fixed-size matrices it could be a CwiseNullaryOp<> thus making the assignment more efficient and optimizable by the compiler. > Is this approach (using a custom scalar type) still the recommended way to create sparse block-matrices? As I said, there is currently no recommended way to create sparse-block matrices but this is something we'd like to have, and we'll probably work on it soon. So if you're willing to share your advances in this direction and participate in any way to this development, you're very welcome!! > Can we reduce the overhead of these initialisations or make them optional? Do you have some measure of the overhead? gael > > I'd welcome any feedback. :) > > Thanks! > Daniel Vollmer > > -------------------------- > Deutsches Zentrum für Luft- und Raumfahrt e.V. (DLR) > German Aerospace Center > Institute of Aerodynamics and Flow Technology | Lilienthalplatz 7 | 38108 Braunschweig | Germany

**Follow-Ups**:**Re: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?***From:*Sameer Agarwal

**RE: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?***From:*Daniel.Vollmer

**RE: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?***From:*Daniel.Vollmer

**References**:

**Messages sorted by:**[ date | thread ]- Prev by Date:
**Re: [eigen] Bug in language_support.cmake** - Next by Date:
**Re: [eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?** - Previous by thread:
**[eigen] SparseMatrix with non POD-"Scalar" type (e.g. SparseBlockMatrix) yields unnecessary initialisations of Scalar-type?** - Next by thread:

Mail converted by MHonArc 2.6.19+ | http://listengine.tuxfamily.org/ |