[petsc-dev] Suggestions for MatProductCreate()

[Lisandro Dalcin]

On Tue, 24 Mar 2020 at 05:02, hzhang at mcs.anl.gov <hzhang at mcs.anl.gov> wrote:

> Matt and Lisandro,
> Let me explain a little more about new API for mat-mat operations, which
> include MatMaMult (AB), MatPtAP (PtAP), MatMatMatMult (ABC), with a total 6
> matproduct types.
>
> 1) MatProduct is logically similar to MatFactor, not a single Mat. Calling
> MatProductCreate(A,B,NULL,&D);
> MatProductSetType(D, MATPRODUCT_AB);
> MatDestory(&A); // A is gone
> MatDestory(&B); // B is gone
> MatProductSymbolic(D);
>  is equivalent to calling MatMatMult(null, null, ..., &D)
> or MatGetFactor(null, ...,&F).
> If you do this, you'll get segfault (I'll add additional error message for
> this use case).
>
>
No, it is not really equivalent to MatMatMult(null, null). Because if you
pass NULL to MatMatMult you can check for them. However, in the new
MatProduct API, if you desdtroy A,B and next call MatProductSymbolic(D) as
in my example, the internal pointers to A and B are not NULL, they are
"dangling" pointers, holding the same addresses at MatProductCreate() time.

This is way I said that this API is prone to dangling pointer issues if
ownership is not properly handled.


All I'm asking is you to implement things this way:

static PetscErrorCode MatProductCreate_Private(Mat A,Mat B,Mat C,Mat D)
{
  ...
  PetscObjectReference(A); // take ownership
  PetscObjectReference(B); // take ownership
  PetscObjectReference(C); // take ownership
  product->A        = A;
  product->B        = B;
  product->C        = C;
}

PetscErrorCode MatProductReplaceMats(Mat A,Mat B,Mat C,Mat D)
{
  ...
  if (A) {
    if (!product->Areplaced) {
      PetscObjectReference(A); // take ownership of input first
      MatDestroy(&product->A); // release old reference
      product->A = A;
    } else ...
  }

  ...
}

PetscErrorCode MatDestroy(Mat *A)
{
 ...
    if ((*A)->product) {
    Mat_Product  *product = (*A)->product;
    MatDestroy(&product->A); // release reference
    MatDestroy(&product->B); // release reference
    MatDestroy(&product->C); // release reference
    MatDestroy(&product->Dwork);
    PetscFree(product);
  }
...
}


> 2) New API is an incremental improvement over the old API, e.g.,
> MatMatMultSymbolic(A,B, fill,&C); MatMatMultNumeric(A,B,C)
> are replaced by
> MatProductCreate(A,B,NULL,&D);
> MatProductSetType(D, MATPRODUCT_AB);
>  MatProductSetFill(D,fill);
>  MatProductSetFromOptions(D);
> MatProductSymbolic(D):
> MatProductNumeric(D);
>
> There are at least two reasons we decided to introduce new API.
> a) there are multiple implementations (algorithms) for some product types
> and matrix types, e.g., A*B  and PtAP for mpiaij_mpiaij matrices.
> Developers are adding more and more special implementations, e.g, for gpu
> etc. Previously, the options of algorithms are selected in
> MatMatMultSymbolic_xxx(...,&D) with random option names. However, some
> options need to be selected before D is setup, and attached to D. For this
> reason, we introduced new API, and add MatProductSetFromOptions(D) to
> handle all the options. Note: the options can only be determined at the
> level of producttype and mattype are known, thus the options are enabled in
> MatProductSetFromOptions_mattype_producttype().
>

As I said in previous emails, I'm not complaining about the API, I said it
was a welcome adition. All I'm asking is:

1)  Fix the handling of reference counting. This is an implementation
detail, not really an API change.

2) Add a simple new routine to the API, MatProductReset(D)  [or maybe name
it MatProductClear(D) ], that should throw away all of the D->product stuff
to release memory related to the product algorithm implementation. After
calling MatProductReset(D), you cannot call anymore MatProductNumeric().



> b) Reuse of MatProduct D = A*B means same A and B with updated values.
> However, some users change A or B, some use different  D without knowing
> original D might contain internal data structure that is associated with A
> and B. The new API forces reuse same A, B, unless user
> calls MatProductReplaceMats(A,B,C,D).
>

All that is good, no objections. But please consider
adding MatProductReset(D) with the behavior I described above.
Of course, this is only needed if you implement the changes to reference
counting that I requested at the beginning.


>
> 3) For MatProduct, I consider
> MatProductSetFromOptions() is logically equivalent to MatSetFromOptions(),
> MatProductSymbolic() is logically equivalent to MatSetUp(),
> and
> MatProductNumeric() assembles the matproduct  D repeatedly when A or B are
> updated.
>
>
OK, so... your point it that what looks like "setup" code
in MatProductSetFromOptions() is just some minimal stuff required to
actually dispatch the SetFromOptions_typeA_typeB routine, right? In that
case, you are right, no changes/refactors needed here.
-- 
Lisandro Dalcin
============
Research Scientist
Extreme Computing Research Center (ECRC)
King Abdullah University of Science and Technology (KAUST)
http://ecrc.kaust.edu.sa/
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20200324/128845d1/attachment-0001.html>