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 _______________________________________________ DOLFIN-dev mailing list [email protected] http://www.fenics.org/mailman/listinfo/dolfin-dev
