On 11 Mar 2019, at 23:59, Jed Brown via petsc-dev <[email protected]<mailto:[email protected]>> wrote:
One could imagine PETSc written so that all Mat operations like transpose and products return lazy representations, with a method that makes it explicit. Would eliminate all the Transpose variants in the matrix products, for example, and enable an implementation to choose order of operations for chained products. That would be cool! This would solve e.g. my concerns related to transpose in MatLoad_HDF5 with PETSC_VIEWER_HDF5_MAT and MATLAB engine interface, see https://bitbucket.org/petsc/petsc/pull-requests/1381/matload_hdf5-w-petsc_viewer_hdf5_mat I was thinking about introducing a special MatOption that would say the (non-preallocated) matrix can be replaced by an implicit one. This wouldn't change the API. But what you suggest sounds much more elegant. Thanks, Vaclav "Smith, Barry F. via petsc-dev" <[email protected]<mailto:[email protected]>> writes: On Mar 11, 2019, at 1:28 PM, Isaac, Tobin G <[email protected]<mailto:[email protected]>> wrote: Okay, but do you have an opinion on the original issue: should we always return a matrix from various matrix-matrix routines, even if it has limited functionality? I'm fine with it always returning the composite matrix unless you have a specific use case that would require MAT_INITIAL_MATRIX_GET_VALUES Barry On Mon, Mar 11, 2019 at 06:09:55PM +0000, Smith, Barry F. wrote: Sounds like overkill on the complexity meter. Barry On Mar 11, 2019, at 6:15 AM, Isaac, Tobin G via petsc-dev <[email protected]<mailto:[email protected]>> wrote: One thought: we could extend MatReuse with MAT_INITIAL_MATRIX_GET_VALUES to indicate that a matrix that implements MatGetValues() is desired, to allow some user intent. On March 9, 2019 7:44:06 PM EST, Matthew Knepley <[email protected]<mailto:[email protected]>> wrote: On Sat, Mar 9, 2019 at 4:42 PM Isaac, Tobin G via petsc-dev < [email protected]<mailto:[email protected]>> wrote: Let's say I want to do Newton-Krylov on some DAG of computations. Building the Jacobian is programmatic from the chain rule, requiring just a lot of matrix-matrix multiplications in the correct sequence. If it really is more efficient for two Jacobians to be contracted into a matrix with random access, that's great, but if it's not, or if the appropriate implementation of MatMatMult_X_Y does not exist, having a MATCOMPOSITE is better than nothing. So why not make MATCOMPOSITE something that our matrix-matrix multiplication routines return when a specialization isn't found? It moves the "Not Implemented" warning from MatMatMult() to whenever I try to do something other than multiply with the output, but in some situations that's desirable. I'd be interested to hear what people think. I think this is a good idea. The most important thing is to maintain transparency, so that the output makes this clear in -ksp_view and we are easily able to see what someone has produced. Matt Thanks, Toby
