[petsc-dev] DMGetMatrix --> DMGetMatrices?

[Dmitry Karpeev]

On Fri, Feb 10, 2012 at 7:45 PM, Jed Brown <jedbrown at mcs.anl.gov> wrote:

> On Fri, Feb 10, 2012 at 19:33, Dmitry Karpeev <karpeev at mcs.anl.gov> wrote:
>
>> How can there be a single callback, when the two DMs could be of
>> different types?
>> I guess it could be done using the double-dispatch mechanism:
>> DMComputeOperators(DM dmJ, DM dmP, Mat *J, Mat *P) would dispatch based
>> on the types of the two DMs,
>> and, if no such method exists, defaulting internally to calling
>> DMComputeMatrix(DM dm, Mat *A) for both J and P.
>>
>
> Or call the user's function which just assembles into the P matrix.
>
>
>> I still think that DMComputeMatrix() should be different from
>> DMCreateMatrix() -- the latter just preallocates --
>> eliminating the need for DMSetPreallocateOnly().
>>
>
> How many times do we have to explain that DMCreateMatrix() *NEVER* fills
> in useful values? It has two capabilities: preallocate (just specify number
> of nonzeros per row) and also setting the column indices. It does both by
> default because that is better for performance profiling and better to
> catch setting values outside of the allocated sparsity pattern at the first
> incident instead of whenever a row overflows. Under no circumstances does
> it assemble physics and under no circumstances does DMComputeJacobian()
> allocate a matrix.
>

How does this explanation square with the following code from KSPSetUp?

 if (ksp->dmActive && !ksp->setupstage) {
    /* first time in so build matrix and vector data structures using DM */
    if (!ksp->vec_rhs) {ierr =
DMCreateGlobalVector(ksp->dm,&ksp->vec_rhs);CHKERRQ(ierr);}
    if (!ksp->vec_sol) {ierr =
DMCreateGlobalVector(ksp->dm,&ksp->vec_sol);CHKERRQ(ierr);}
    ierr = DMCreateMatrix(ksp->dm,MATAIJ,&A);CHKERRQ(ierr);
    ierr = KSPSetOperators(ksp,A,A,stflg);CHKERRQ(ierr);
    ierr = PetscObjectDereference((PetscObject)A);CHKERRQ(ierr);
 }

Dmitry.


>
>>
>> If the two DMs are of the same type and, say, share a mesh, then the dual
>> dispatch can take advantage of
>> a corresponding implementation.
>>
>
> The user should always be able to get the callback to do both.
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20120210/538acf8e/attachment.html>