> On Jul 29, 2017, at 7:20 AM, Karl Rupp <r...@iue.tuwien.ac.at> wrote:
> 
> 
>>>>>       In theory PETSc Mat have two "state" values,
>>>>> 
>>>>> 1) nonzerostate - this is increased anytime the nonzero structure changes
>>>>> 
>>>>> 2) state - this is increased anytime any numerical values are changed
>>>>> 
>>>>>    These are used by PCSetUp() to determine if the preconditioner needs 
>>>>> to be updated and if that involves a new nonzero structure.
>>>>> The reason I say "in theory" is that it looks like certain changes to 
>>>>> matrices do not properly update the state.
>>>>> 
>>>>> Hi Barry: Sounds good. In the meantime, I'm simply making a corresponding 
>>>>> _SeqAIJMKL version of all of the functions in question that will just 
>>>>> call the _SeqAIJ version and then call a function that creates an updated 
>>>>> MKL sparse matrix handle. It appears that several of these do not have 
>>>>> function prototypes in any of the .h files (they haven't been needed 
>>>>> outside of aij.c). I assume I can just add PETSC_INTERN prototypes in 
>>>>> src/mat/impls/aij/seq/aij.h (where we have things like prototypes for 
>>>>> MatAssemblyEnd_SeqAIJ)?
>>>>   I think this is the correct place. If they have static in front of them 
>>>> you will need to remove that also.
>>> 
>>> I'm just wondering: Isn't fixing the state variables the faster, easier and 
>>> more maintainable approach?
>>  So you fix the state variable by reseting with each each function that 
>> changes the matrix entries, how does that magically cause the MKL "convert 
>> to MKL format" in these cases? It is a slightly orthogonal issue I think.
> 
> 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?

   For parallel matrices are all the operations that need to do this kind of 
update collective over the Mat?

> 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.

  Barry


> 
> 
>>   Or for every PETSc object do we optionally provide a callback function 
>> that is called with each change to state? Currently this would not be 
>> efficient because there may be several (many) changes to state before we 
>> want the subtype data structure to be updated. At the moment we have a 
>> "hard" MatAssemblyEnd() that handles changes in nonzero structure, may be we 
>> need something similar for changes in numerical value? BTW this is also true 
>> for CUDA/CL representations?
> 
> A callback for every state change is too greedy, yes. As long as state and 
> nonzerostate are reliable, we should be good in terms of outside dependencies 
> on the matrix values and structure, respectively.


> 
> Best regards,
> Karli

Reply via email to