[petsc-users] Extracting data from a Petsc matrix

[Barry Smith]

On Jul 16, 2013, at 2:18 PM, Harshad Sahasrabudhe <hsahasra at purdue.edu> wrote:

> Hi Barry,
> 
> I'm confused, can you please explain what you meant by 'MAGMA has a calling sequence like LAPACK'?

Here is how we call lapack Cholesky factorization from PETSc:

#undef __FUNCT__
#define __FUNCT__ "MatCholeskyFactor_SeqDense"
PetscErrorCode MatCholeskyFactor_SeqDense(Mat A,IS perm,const MatFactorInfo *factinfo)
{
#if defined(PETSC_MISSING_LAPACK_POTRF)
  PetscFunctionBegin;
  SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"POTRF - Lapack routine is unavailable.");
#else
  Mat_SeqDense   *mat = (Mat_SeqDense*)A->data;
  PetscErrorCode ierr;
  PetscBLASInt   info,n;

  PetscFunctionBegin;
  ierr = PetscBLASIntCast(A->cmap->n,&n);CHKERRQ(ierr);
  ierr = PetscFree(mat->pivots);CHKERRQ(ierr);

  if (!A->rmap->n || !A->cmap->n) PetscFunctionReturn(0);
  PetscStackCallBLAS("LAPACKpotrf",LAPACKpotrf_("L",&n,mat->v,&mat->lda,&info));
  if (info) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_MAT_CH_ZRPVT,"Bad factorization: zero pivot in row %D",(PetscInt)info-1);
  A->ops->solve             = MatSolve_SeqDense;
  A->ops->solvetranspose    = MatSolveTranspose_SeqDense;
  A->ops->solveadd          = MatSolveAdd_SeqDense;
  A->ops->solvetransposeadd = MatSolveTransposeAdd_SeqDense;
  A->factortype             = MAT_FACTOR_CHOLESKY;

  ierr = PetscLogFlops((A->cmap->n*A->cmap->n*A->cmap->n)/3.0);CHKERRQ(ierr);
#endif
  PetscFunctionReturn(0);
}


   If Magma has a calling sequence similar to LAPACKpotrf_("L",&n,mat->v,&mat->lda,&info)); then you could model your code on the routine above instead of typing in all new code that does almost the exact same thing.

   Barry

> 
> Thanks,
> Harshad
> 
> On 07/16/2013 02:48 PM, Barry Smith wrote:
>>    Read all of http://www.mcs.anl.gov/petsc/developers/index.html
>> 
>>     Note that if Magma has a calling sequence like lapack you could possible steal chunks of code from the routines I pointed you to yesterday and modify them as needed so you don't need to reinvent the wheel.
>> 
>> 
>>    Barry
>> 
>> On Jul 16, 2013, at 1:37 PM, Matthew Knepley <knepley at gmail.com> wrote:
>> 
>>> On Tue, Jul 16, 2013 at 1:13 PM, Harshad Sahasrabudhe <hsahasra at purdue.edu> wrote:
>>> Hi Jed,
>>> 
>>> Thanks for your reply.
>>> 
>>> You're on your own for storage of factors.  Alternatively, you could add
>>> library support so that you could use PCLU and
>>> '-pc_factor_mat_solver_package magma' (or PCFactorSetMatSolverPackage).
>>> Doing this is not a priority for us, but we can provide guidance if you
>>> want to tackle it.
>>> 
>>> I would definitely like to start working on adding library support. I think this is the most efficient way to go about it. Can you give me certain details such as:
>>> 
>>> 1) How should I start going about it?
>>> 
>>> Read the UMFPACK implementation
>>>  2) How will I check-in the changes to Petsc?
>>> 
>>> Using Git
>>>  3) What version of Petsc will the changes be reflected in if I started working on it right now?
>>> 
>>> A branch of 'master'
>>>  4) How many hours does it generally take to get this done?
>>> 
>>> How many licks does it take to get to the center of a Tootsie Roll Pop?
>>>  5) How is the peer review done?
>>> 
>>> Through a pull request on BitBucket.
>>> 
>>>   Thanks,
>>> 
>>>      Matt
>>>  Thanks,
>>> Harshad
>>> 
>>> On 07/13/2013 12:43 PM, Jed Brown wrote:
>>> "hsahasra at purdue.edu" <hsahasra at purdue.edu> writes:
>>> 
>>> Hi,
>>> 
>>> I am working on solving a system of linear equations with square
>>> matrix. I'm first factoring the matrix using LU decomposition.
>>> I assume you're solving a dense problem because that is all MAGMA does.
>>> 
>>> I want to do the LU decomposition step using MAGMA on GPUs. MAGMA
>>> library implements LAPACK functions on a CPU+GPU based system.
>>> 
>>> So my question is, how do I extract the data from a Petsc Mat so that
>>> it can be sent to the dgetrf routine in MAGMA.
>>> MatDenseGetArray
>>> 
>>> Is there any need for duplicating the data for this step?
>>> You're on your own for storage of factors.  Alternatively, you could add
>>> library support so that you could use PCLU and
>>> '-pc_factor_mat_solver_package magma' (or PCFactorSetMatSolverPackage).
>>> Doing this is not a priority for us, but we can provide guidance if you
>>> want to tackle it.
>>> 
>>> 
>>> 
>>> 
>>> -- 
>>> What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.
>>> -- Norbert Wiener
>