| Thread Previous • Date Previous • Date Next • Thread Next |
Johan Hake wrote:
On Thursday 05 March 2009 20:36:06 Garth N. Wells wrote:Anders Logg wrote:On Thu, Mar 05, 2009 at 08:12:29PM +0100, Johan Hake wrote:On Thursday 05 March 2009 17:13:02 Garth N. Wells wrote:We have a high-level interface which provides operations like A += B for matrices.For the record: This is not only for the python interface. The c++ operator+= also use axpy.This causes some problems as it hides the complexity which is inherent in sparse matrices, in particularly with respect to whether or not the two matrices have the same non-zero pattern. Unless someone has a good idea as to get around this robustly, I suggest that we remove these functions for matrices. A user can still do A.axpy(1.0, B) or A.axpy(1.0, B, True) where the last argument indicates whether or not the sparsity patterns of A and B are the same (default is False).A naive suggestion: Is it possible to compute a unique number based on our SparsityPattern? This number could be stored as a private member of a Matrix? Then only matrices from the same SparsityPattern will be addable. When ever a Matrix is changed, by other means than M.init(sp), this number will be set to some default incompatible number. If we do not find a robust way to do this I am fine with removing the operators that use axpy, but I really think it is a neat feature which I would like to keep.How about having both axpy() and operators like +, += but implemented differently. The axpy() call can make assumptions about the sparsity pattern (or take a flag) and use library calls to PETSc and Epetra to compute the result.I would hope that we can make axpy() do the work for +=, =+, etc, with the assumption that the sparsity patterns do differ. The tricky bit is getting this to work with some of the backends, namely Epetera which appears to me to be less well documented than PETSc so I couldn't figure it out easily.Have you experienced that two matrices that do not have the same nonzero pattern are not addable using the EpetraExt::MatrixMatrix::Add? At least nothing of that sort is mentioned in the Doxygen generated documentation.
It's mentioned here http://trilinos.sandia.gov/packages/docs/dev/packages/epetraext/doc/html/classEpetraExt_1_1MatrixMatrix.html#e1
If this is right we only have the "problem" for the PETSc backend. Then we can change back the GenericMatrix interface to an axpy without same_nonzero_patter. The axpy function in PETScMatrix can use MatAXPY with DIFFERENT_NONZERO_PATTERN as default, and then we can add a second axpy taking the same_nonzero_pattern as argument.
Let's keep it as is because even for backends where 'same_nonzero_pattern' its not yet used, it could be used the future to get some extra performance.
Garth
Johan
| Thread Previous • Date Previous • Date Next • Thread Next |
