[petsc-users] A number of questions about DMDA with SNES and Quasi-Newton methods

[Jed Brown]

Hmm, this is a use case we would like to support, but I'm not sure how
to do it.  I think at this point that your local submatrix will only
have a stencil set if you were using MatNest, which I don't recommend
(whatever you do, don't depend on it).  But since you have the DMDA for
your J_hh block, you can call MatSetStencil() yourself.

  ierr = DMDAGetGhostCorners(da,&starts[0],&starts[1],&starts[2],&dims[0],&dims[1],&dims[2]);CHKERRQ(ierr);
  ierr = MatSetStencil(A,dim,dims,starts,dof);CHKERRQ(ierr);

I don't know of code you could use to translate the stencil coordinates
to local indices (it's just a lexicographic ordering so simple for you
to write), but you can use MatSetValuesLocal on the J_hr/J_rh blocks.

"zakaryah ." <zakaryah at gmail.com> writes:

> Well I made a little progress by considering SNES ex28.c.  In the Jacobian
> routine, I call DMCompositeGetLocalISs, then use the IS to call
> MatGetLocalSubMatrix.  I call these J_rr, J_rh, J_hr, and J_hh, where r
> represents the redundant variables and h represents the displacements.  I
> assume I can call MatSetValuesStencil on J_hh, as before, and MatSetValue
> on J_rr (which is 1x1).  I'm guessing that J_rr, J_rh, and J_hr can only be
> set on the processor which owns the redundant variable - is this correct?
> How do I determine the ordering for J_hr and J_rh?
>
> On Tue, Oct 24, 2017 at 2:45 PM, zakaryah . <zakaryah at gmail.com> wrote:
>
>> I see - I use a local variable to add up the terms on each process, then
>> call MPI_Reduce within the function on process 0, which owns the redundant
>> variable.
>>
>> I have one more question - for the calculation of the Jacobian, my life is
>> made much much easier by using MatSetValuesStencil.  However, the matrix
>> which the SNES Jacobian receives as argument is the "full" matrix,
>> containing the DMDA variables (displacements), plus the redundant
>> variable.  How do I access the submatrix corresponding just to the DMDA?
>> If I can do that, then I can call MatSetValuesStencil on the submatrix.  Is
>> this the right approach?  I'm not sure how to set the elements of the
>> Jacobian which correspond to the redundant variable, either - i.e., how do
>> I get the ordering?
>>
>>