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. The other calls can use getrow() / setrow() and do everything through the GenericMatrix interface (so no separate implementation is needed for each backend). This will be slower than axpy() but it is very convenient to write A + B and sometimes speed is not important. If it is, then use axpy(). -- Anders
signature.asc
Description: Digital signature
_______________________________________________ DOLFIN-dev mailing list [email protected] http://www.fenics.org/mailman/listinfo/dolfin-dev
