On Sep 18, 2013, at 10:52 AM, "Vijay S. Mahadevan" <[email protected]> wrote:
>> Why are you setting this flag? If you don't set it then the code will go
>> through.
>
> Yes it passes through with DMDA. But with my custom DM implementation
> (MOAB based), this fails unless I explicitly set the NNZ pattern in
> DMCreateMatrix. I currently don't do that and the only way to
> replicate this behavior with DMDA was setting the flag explicitly.
Until you can provide the preallocation you need to simply call MatSetUp()
and then MatAssemblyBegin/End() on your matrix inside your
DMCreateMatrix_Moab(). The expectation is that the DM matrix provided is "ready
to go", we didn't expect to handle matrices that were "incomplete" in other
parts of the code.
>
>> There is no reason to use a user created matrix since the DM can provide
>> the correct one. One of the "main jobs" of DM is to create that matrix so we
>> really don't intend people to use a DM but then create the matrix themselves.
>
> Agreed, and that's what I expected when I passed in NULL args instead
> of user matrix. But the DM created matrix errors out during Duplicate
> while user provided matrix passes through, even though both of them
> are unassembled.
PetscErrorCode MatDuplicate(Mat mat,MatDuplicateOption op,Mat *M)
{
PetscErrorCode ierr;
Mat B;
PetscInt i;
PetscFunctionBegin;
PetscValidHeaderSpecific(mat,MAT_CLASSID,1);
PetscValidType(mat,1);
PetscValidPointer(M,3);
if (!mat->assembled)
SETERRQ(PetscObjectComm((PetscObject)mat),PETSC_ERR_ARG_WRONGSTATE,"Not for
unassembled matrix");
The "user matrix" must have been assembled otherwise it would have stopped
above.
> And this is what is confusing.
>
> Vijay
>
> On Wed, Sep 18, 2013 at 10:46 AM, Barry Smith <[email protected]> wrote:
>>
>> The code did not "crash", it ended with an error message
>>
>>
>> On Sep 18, 2013, at 10:32 AM, "Vijay S. Mahadevan" <[email protected]> wrote:
>>
>>> All,
>>>
>>> When I call DMSetMatrixPreallocateOnly
>>
>> Why are you setting this flag? If you don't set it then the code will go
>> through.
>>
>>> and then TSSetRHSJacobian
>>> without a matrix input argument (NULL instead), then the code crashes
>>> during TSSolve and more specifically during MatDuplicate. This is
>>> replicable with src/ts/examples/tutorials/ex2.c. Here's the error with
>>
>> Jed is reworking the logic here so the code will change.
>>
>>
>>> stack trace.
>>>
>>> [0]PETSC ERROR: --------------------- Error Message
>>>
>>>
>>> I hope this is not supposed to happen. The fix seems to be to call
>>> TSSetRHSJacobian with a user created matrix instead of letting DM
>>> create one,
>>
>> There is no reason to use a user created matrix since the DM can provide
>> the correct one. One of the "main jobs" of DM is to create that matrix so we
>> really don't intend people to use a DM but then create the matrix themselves.
>>
>> Barry
>>
>>
>>> which rectifies this issue. This behavior feels
>>> inconsistent and is there a way to fix this.
>>>
>>> Vijay
>>>
>>> PS: Just replace the following lines in ex2.c to replicate the error.
>>>
>>> /* ierr = TSSetRHSJacobian(ts,A,A,RHSJacobian,&appctx);CHKERRQ(ierr); */
>>> ierr = TSSetRHSJacobian(ts,NULL,NULL,RHSJacobian,&appctx);CHKERRQ(ierr);
>>> ierr = DMSetMatrixPreallocateOnly(appctx.da,PETSC_TRUE);CHKERRQ(ierr);
>>