[eigen] Re: Special Matrices: first get diagonal matrices done |

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

*To*: eigen <eigen@xxxxxxxxxxxxxxxxxxx>*Subject*: [eigen] Re: Special Matrices: first get diagonal matrices done*From*: Benoit Jacob <jacob.benoit.1@xxxxxxxxx>*Date*: Sun, 28 Jun 2009 21:56:35 +0200*Dkim-signature*: v=1; a=rsa-sha256; c=relaxed/relaxed; d=gmail.com; s=gamma; h=domainkey-signature:mime-version:received:in-reply-to:references :date:message-id:subject:from:to:content-type :content-transfer-encoding; bh=ZNFi4wpoJ7xCSM+Q2qAV62SxTl3nEWB7zu+zKul6Yw4=; b=XZLXNR33FAj1QMFhS7iz5YI+k0djrIHstfjtwAe5DQVuIRBFzbEE/4G3F0Vz3m9loE fEGCUuupe0YKIDN9GpTssYJdtEBIjK4G7yFSE2OHILVa48D2MuK/OTDavkPXa0PZP1Jj WbIpNtJE69G4ZliCN6Q8s7b7gpTculTQgP0Oo=*Domainkey-signature*: a=rsa-sha1; c=nofws; d=gmail.com; s=gamma; h=mime-version:in-reply-to:references:date:message-id:subject:from:to :content-type:content-transfer-encoding; b=PvGOhTQwDPLwkpKU00lmTGc/UlZ2IpJ7J/7vUkq3jE+xgIwQBs6n/BZBX9wTPbBR8H fD8JoRUCl4Jb9avJxrNu5fDRgnZIuX8/uIFWg6XRfk0PGRJR7uxoNDcZqWfur9i0pOBG 1nvjnl6QKaw6xzHCNG07bYPbc71ByjN4paZMk=

Hi, Finally, this is done, in this fork: http://bitbucket.org/bjacob/eigen2_special_matrices/ I imported permissions from eigen/eigen2 so you can commit directly to it. The main reasons why I don't push directly to eigen/eigen2 are that it's quite invasive, you may have comments, and also I was too lazy to update the Sparse module accordingly, so I hope that Gael can help me here; what needs to be done seems mainly to update SparseDiagonalProduct.h, in the same way as the file Core/DiagonalProduct.h: decouple it completely from the other products by implementing a separate operator* for DiagonalBase. So let's explain the new class hierarchy. I did some renaming: DiagonalMatrixBase--->DiagonalBase DiagonalMatrixWrapper--->DiagonalWrapper DiagonalBase is still the common parent of DiagonalMatrix and of DiagonalWrapper. What's new is 1) DiagonalBase doesn't inherit MatrixBase anymore, and 2) there are new "abstract" base classes encompassing both DiagonalBase and MatrixBase. These new classes are MultiplierBase and AnyMatrixBase. See file MatrixBase.h. A AnyMatrixBase is anything that can be copied into a MatrixBase. A MultiplierBase is a AnyMatrixBase that furthermore can be multiplied with a MatrixBase, producing a new MatrixBase. AnyMatrixBase and MultiplierBase just serve to overload operator= , operator*, etc. The only method in them is derived(). DiagonalBase() is very minimalistic, for example I didn't implement rows(), cols() etc. Instead, it has a method diagonal() that gives the actual vector of diagonal coefficients. The idea was that at best we would implement a very limited subset of the MatrixBase API in DiagonalBase. That would force the user to look up the docs frequently for availability of methods in DiagonalBase. So instead, I thought that it would be easiest to force the user to address explicitly the vector of diagonal coeffs, diagonal(), and from that point on the user is in known territory with an actual MatrixBase. That said we can always add rows() etc if you think that's useful. The next steps for me are to code BandedMatrix (should now be a piece of cake) and start the SVD... Gael--->if you can take care of adapting the Sparse module, that would be wonderful on another note, I've been thinking too that the whole marked() stuff for triangular and selfadjoint matrices could go away, in fact I don't believe anymore that triangularness should be part of the Flags. Instead it would be simpler and IMO just as good to have separate methods like selfAdjointProduct(). Cheers, Benoit 2009/5/11 Benoit Jacob <jacob.benoit.1@xxxxxxxxx>: > Hi, > > Before getting started with BandedMatrix and other kinds of special > matrices, I think the perfect starting point is to first nail the > situation with diagonal matrices. > > Currently we have: > * class Diagonal nests a matrix xpr and is a vector xpr > * class DiagonalMatrixWrapper nests a vector xpr and is a matrix xpr > * class DiagonalMatrix stores a vector and is a matrix xpr > > I think that we agree that we want instead DiagonalMatrix to NOT be a > matrix xpr (ie NOT inherit MatrixBase) and have a coeff() method that > just asserts row==col. > > DiagonalMatrixWrapper is still as much needed as before for situations like > > mat * vec.asDiagonal(); > > where asDiagonal() is actually just being used as a thin wrapper > around the vector allowing to select the matrix product strategy. > > What is open for discussion is whether DiagonalMatrixWrapper should > still inherit MatrixBase? And if not, should there be an AsDense kind > of class constructing a MatrixBase expression out of a > DiagonalMatrixWrapper without evaluating? > > Let me explain why I think that the answer to both questions is _no_. > The solution i'm going to advocate will be to keep the same class > names, only make DiagonalMatrixBase no longer inherit MatrixBase. > > First of all we could have a DiagonalMatrixBase::toMatrix() that > returns a plain matrix (see at the end for naming discussion), so that > the user who wants something that behaves like a matrix, can get a > plain matrix if he wants to: > > cout << vec.asDiagonal().toMatrix() << endl; > > The argument for having instead a asDense() returning a MatrixBase, is > to avoid the cost of evaluating the matrix, but then remember that > this expression would have a very expensive coeff() method as it would > have to do a if(row==col). So as soon as you're accessing a large part > of the matrix, it's not good anyway (performance wise) to have an > asDense expression here. > > I also think that there isn't much of a use case for doing directly > "vec.asDiagonal()(i,j)" -- it's OK, i think, to force the user to > replace that by "i==j?vec(i):0" making him realize where there is room > for optimization. > > Then let's look at the main use case for asDiagonal(): > > mat * vec.asDiagonal() > > Currently, since DiagonalMatrixWrapper is a MatrixBase, we have to > catch that in Product.h, adding complexity in there. > > If instead DiagonalMatrixWrapper were NOT a MatrixBase, we could catch > that in a completely separate operator*. Less code overall and > especially in Product.h, shorter compilations times, "you don't pay > for what you don't use", same result. > > Is there anything I didn't think of? > > Yes that means that one could no longer do: > mat * mat2.marked<DiagonalBits>() > and hope for that product to be optimized, but frankly !! if you know > that a matrix happens to be diagonal, then you can probably avoid > storing it as a full matrix in the first place! So I don't think > there's a real use case for that; and even if there were one could > still do: > mat * mat2.diagonal().asDiagonal(). > > To summarize my plan: > *no more changes to Diagonal (the ex-DiagonalCoeffs) > *make DiagonalMatrixBase no longer inherit MatrixBase. It doesn't even > need a coeff() method. > *DiagonalMatrix still inherits DiagonalMatrixBase > *DiagonalMatrixWrapper still inherits DiagonalMatrixBase > *DiagonalMatrixBase is the class against which we implement the > operator* catching diagonal products, and by the way we can also have > operator= etc. > *DiagonalMatrixBase::toMatrix() returns a plain Matrix. > > (the last method toMatrix(), i couldn't find a nicer name: eval() > suggests it returns a DiagonalMatrix, toDense() doesn't make it clear > it evaluates, ...) > > Notice that if you are really in favor of a diagonal matrix > expression, we can still add this to this scheme using multiple > inheritance on MatrixBase and DiagonalMatrixBase... > > Cheers, > Benoit >

**Messages sorted by:**[ date | thread ]- Prev by Date:
**Re: [eigen] Std Vector again** - Next by Date:
**[eigen] Non linear equations solver** - Previous by thread:
**Re: [eigen] quaternion transcendentals** - Next by thread:
**[eigen] Non linear equations solver**

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