Well, I thought of AIJMKL keeping track of the state for which
mkl_sparse_optimize() was called. If the matrix state changes, the next call to
MatMult()
is MatMult the only operation that needs this check, or do many Mat methods
need this check?
Looks like it is needed for MatMult(), MatMatMul() and triangular solves:
https://software.intel.com/en-us/mkl-developer-reference-fortran-inspector-executor-sparse-blas-analysis-routines
For parallel matrices are all the operations that need to do this kind of
update collective over the Mat?
According to the list above, all operations are collective.
will detect that the current state does not match the state it was optimized
for and hence trigger another optimization. This isn't too different from what
we do for GPU stuff.
We do this for cusparse
static PetscErrorCode MatMult_SeqAIJCUSPARSE(Mat A,Vec xx,Vec yy)
{
....
PetscFunctionBegin;
/* The line below is necessary due to the operations that modify the matrix
on the CPU (axpy, scale, etc) */
ierr = MatSeqAIJCUSPARSECopyToGPU(A);CHKERRQ(ierr);
But for CUSP (deprecated) we do it at MatAssembly time
PetscErrorCode MatAssemblyEnd_SeqAIJCUSP(Mat A,MatAssemblyType mode)
{
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = MatAssemblyEnd_SeqAIJ(A,mode);CHKERRQ(ierr);
ierr = MatCUSPCopyToGPU(A);CHKERRQ(ierr);
as we do it for ViennaCL
PetscErrorCode MatAssemblyEnd_SeqAIJViennaCL(Mat A,MatAssemblyType mode)
{
PetscErrorCode ierr;
PetscFunctionBegin;
ierr = MatAssemblyEnd_SeqAIJ(A,mode);CHKERRQ(ierr);
ierr = MatViennaCLCopyToGPU(A);CHKERRQ(ierr);
For SeqAIJPERM we do that extra processing of the data in MatAssemblyEnd and
don't need to do anything if only the numerical values change.
So we seem to have some inconsistencies in how we handle this which will
lead to more errors.
ViennaCL should adopt the CUSPARSE model then (and CUSP be removed).
Best regards,
Karli