[petsc-dev] trying to make MKL work with GAMG

Richard Tran Mills rtmills at anl.gov
Fri Jul 6 19:00:18 CDT 2018


On Fri, Jul 6, 2018 at 2:26 PM, Mark Adams <mfadams at lbl.gov> wrote:

>
>> Glad to hear that you got something (at least sort of) working, Mark. Let
>> me know if you want me to work with your users on understanding the
>> performance they are seeing. When I was working for Intel, I probably spent
>> more time on dealing with stupid thread affinity settings than anything
>> else when "optimizing" code performance.
>>
>>
> That is a good idea. I will loop you into the thread.
>
>
>> I should point out that, even with all of the thread settings correct,
>> the current PETSc code base might not get much out of OpenMP threading with
>> MKL. I'm assuming since you are using GAMG that a significant portion of
>> time is spent in MatPtAP() operations? I'm not sure how much use of the
>> sequential matrix-matrix multiply code the parallel MatPtAP() code can
>> achieve.
>>
>
> The RAP or AP is in MKL right? I understand they have these sparse
> matrix-matrix kernels. The get some reuse out of the RAP I think. I will
> ask when I bring you in the thread (but I will wait for you to clarify what
> is supported).
>

My SeqAIJMKL routines currently implement MatMatMult(),
MatMatMultNumeric(), MatTransposeMatMult(), MatPtAP(), and
MatPtAPNumeric(). Notice that there are no MatMatMultSymbolic and
MatPtAPSymbolic() routines in there. The problem is that, even though MKL
internally computes a symbolic multiplication and then a numeric one, until
version 18 update 2 (the latest MKL), there was only one interface
provided: one that did the full (symbolic + numeric) multiplication, that
is, that maps to MatMatMult(). I (and I think a few others) complained
about this, so in version 18 update 2, they added a "two stage" interface.
Unfortunately, it's not the two stages that I want: The first stage only
determines the amount of memory needed but doesn't compute the sparsity
pattern. I did add some "Numeric" routines that use the sparsity pattern
from a previous full multiplication (so for the MAT_REUSE_MATRIX case), but
I've recently learned that what I've done isn't actually safe because that
causes a memory leak in MKL. Frustrating.

Support for a MatPtAP() operation was added very recently in MKL, but it is
only for a symmetric matrix A. I check for MatIsSymmetricKnown(), and call
the MKL version only if this check returns true; otherwise, I just call
MatPtAP_SeqAIJ_SeqAIJ().

--Richard


>
>> There are some issues with mapping between the MKL interfaces for this
>> and the PETSc ones, e.g., there is currently no equivalent of
>> MatMultSymbolic() in MKL. I've had some discussions with the MKL team about
>> this, but the current model they have isn't one that I can map all of the
>> PETSc matrix-matrix interface onto. Hopefully I can get them to change this.
>>
>
> OK, I don't think it will be a problem in the near term ... and I don't
> understand why there is no equivalent of MatMultSymbolic (setting up
> communication maps?). How does if function without this. Please clarify!
>
> Thanks,
>
>
>>
>> --Richard
>>
>>
>>>
>>>>
>>>>
>>>>    Barry
>>>>
>>>>
>>>> > On Wed, Jul 4, 2018 at 9:44 AM Mark Adams <mfadams at lbl.gov> wrote:
>>>> >
>>>> >
>>>> > On Wed, Jul 4, 2018 at 3:01 AM Richard Tran Mills <rtmills at anl.gov>
>>>> wrote:
>>>> > Hi Mark,
>>>> >
>>>> > I'd never looked at the code for MatCreateMPIAIJWithSeqAIJ(), but it
>>>> looks like if you are using MPIAIJ matrices but you've set
>>>> "-mat_seqaij_type seqaijmkl",
>>>> >
>>>> > No, I am using -mat_type aijmkl
>>>> >
>>>> > then the "diagonal" portion of the MPIAIJ matrix ends up as a
>>>> SEQAIJMKL instance, but the off-diagonal portion ends up as a plain SEQAIJ
>>>> when it is created with MatCreateSeqAIJWithArrays(). (I have not stepped
>>>> through this with a debugger to verify that this is actually what happens.)
>>>> So, if we want the off-diagonal matrix to also be a SEQAIJMKL, then we
>>>> either need to add some logic after the MatCreateSeqAIJWithArrays() call to
>>>> check to see if a different subtype of SEQAIJ is being used and then
>>>> convert it to the correct type, or we need a variant of
>>>> MatCreateSeqAIJWithArrays() that will honor the subtype of SEQAIJ that we
>>>> are using. Barry, since you added the concept of subtyping for SEQAIJ, what
>>>> do you prefer?
>>>> >
>>>> > No preference.
>>>> >
>>>> > I note that in some (most?) cases, I think we are OK or even better
>>>> off with the off-diagonal matrix in the MPIAIJ not being a SEQAIJMKL one,
>>>> because the PETSc routines can used the "compressed row" stuff to speed up
>>>> handling the many zero rows that often end up in the off-diagonal, whereas
>>>> the MKL routines cannot.
>>>> >
>>>> >
>>>> > OK , I'll keep that in mind. I pushed a branch with fixes to get my
>>>> gamg tests working (ksp ex56 primarily). I'm going to make a new branch and
>>>> try "-mat_seqaij_type seqaijmkl".
>>>> >
>>>> > Thanks,
>>>> > Mark
>>>> >
>>>> > --Richard
>>>> >
>>>> > On Tue, Jul 3, 2018 at 4:31 PM, Mark Adams <mfadams at lbl.gov> wrote:
>>>> > I have having to fix AIJ methods that don't get the type from the
>>>> args to set created matrix type, so as to keep the MKL typing.
>>>> >
>>>> > I am not sure how to fix this because I don't know to say 'make MPI
>>>> from SEQ' in this method
>>>> >
>>>> > PetscErrorCode MatCreateMPIAIJWithSeqAIJ(MPI_Comm comm,Mat A,Mat
>>>> B,const PetscInt garray[],Mat *mat)
>>>> > {
>>>> >  ....
>>>> >   ierr = MatCreate(comm,mat);CHKERRQ(ierr);
>>>> >   ierr = MatGetSize(A,&m,&n);CHKERRQ(ierr);
>>>> > ....
>>>> >   ierr = MatSetSizes(*mat,m,n,PETSC_DECIDE,N);CHKERRQ(ierr);
>>>> >   ierr = MatSetType(*mat,MATMPIAIJ);CHKERRQ(ierr);
>>>> >
>>>> > How should I do this?  Other places in the code I get the type from
>>>> input (MPI) args and set created matrix types accordingly. I could just put
>>>> a switch with an error for any new types that might come along in the
>>>> future ...
>>>> >
>>>> > Mark
>>>> >
>>>>
>>>>
>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20180706/fd1d4c4f/attachment.html>


More information about the petsc-dev mailing list