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

[Richard Tran Mills]

On Thu, Jul 5, 2018 at 8:32 AM, Mark Adams <mfadams at lbl.gov> wrote:

>
>
> On Wed, Jul 4, 2018 at 1:09 PM Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
>
>>
>>
>> > On Jul 4, 2018, at 9:40 AM, Mark Adams <mfadams at lbl.gov> wrote:
>> >
>> > "-mat_seqaij_type seqaijmkl" just worked.
>>
>>     Mark,
>>
>>       Please clarify. Does this mean you can use -mat_seqaij_type
>> seqaijmkl to satisfy all your needs now without changing any code?
>
>
> Apparently yes. I timed the code and it took much longer with
> NUM_OMP_THREADS > 1 (not sure what that is about but not my problem now),
> so that tells me that threads were activated.
>
> I want to package up this simple approach and my branch that works with -mat_type
> aijmkl , and let the user analyze each. (So I am not worrying about
> performance, I just want correctness).
>

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.

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

--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/3905bf09/attachment-0001.html>