# [petsc-users] Accessing 'halo' matrix entries?

Matthew Knepley knepley at gmail.com
Thu Jun 11 13:41:25 CDT 2015

```On Thu, Jun 11, 2015 at 9:40 AM, Oliver Henrich <ohenrich at epcc.ed.ac.uk>
wrote:

> Dear Matt,
>
> Thanks for your quick response.
>
> The way I understand it is that I can only modify the halo of a vector,
> which is simply copied from the other side of the mesh, but matrices don’t
> actually have halo entries or the like.
>
> I am not sure the solution to this problem is simply adding the difference
> dpsi at one boundary only and not subtracting it at the other side. This
> would mean that at one side you have periodicity plus the added difference,
> whereas at the other side you just have periodicity. Sounds like it leads
> to a kind of discontinuity.
>
> It is also not possible to just add it at the end. It had to be added at
> every step during the iteration process.
>
> Would it be possible to model the halo region explicitly by adding two
> additional points at the leftmost and rightmost boundary, but only on the
> processes with minimal and maximal Cartesian rank along this dimension
> where the jump occurs? Then I could modify the matrix and rhs and impose
> the constraint directly.
>

I think you are missing my point. Perhaps we need an explicit example.
Suppose we have a 1D DMDA,

0 -- 1 -- 2 -- 3

and it is periodic. DMDA by definition has a collocated discretization, so
all unknowns live on the vertices.
If the stencil is size 1, then at vertices 1 we have

0 -- 1 -- 2

which is normal, and at vertex 0 we have

3 -- 0 -- 1

Now this halo 3 should be dp away from 3 itself, since it the periodic
image, so you use p(3) - dp. For the
stencil of 3 we have

2 -- 3 -- 0

Now I have no idea how you make sense of 0 here. You could use p(0) + dp,
or just some Neumann condition here.

I am not convinced this is sensible in higher dimensions since what you
really want to fix is the average pressure
drop which is a much different thing. Thats an integral condition which you
should handle by projection (I think).

Matt

> Best wishes and many thanks for your help,
> Oliver
>
> On 11 Jun 2015, at 14:44, Matthew Knepley <knepley at gmail.com> wrote:
>
> On Thu, Jun 11, 2015 at 8:15 AM, Oliver Henrich <ohenrich at epcc.ed.ac.uk>
> wrote:
>
>> Dear PETSc-Team,
>>
>> I am trying to solve a Poisson equation with a mixed periodic-Dirichlet
>> boundary condition. What I have in mind is e.g. a compressible flow with a
>> total pressure difference imposed between the two sides of the system, but
>> otherwise periodic, and periodic boundary conditions along the remaining
>> two dimensions. Another example would be an electrostatic system with
>> dielectric contrast in an external electric field / potential difference.
>>
>> For clarity, if x = 0  (N+1) is the left (right) halo site at the
>> boundary and x = 1 (N) is the leftmost (rightmost) site in the physical
>> domain:
>>
>> psi(x = 0)  =  psi(x = N) - dpsi
>> psi(x = N+1) = psi(1) + dpsi
>>
>> I know it is possible to solve this with a double Poisson solve, which I
>> try to avoid for performance reasons.
>>
>> It is also possible to solve this by modifying the matrix with a
>> master-slave approach that imposes the constraint. This requires defining a
>> transformation matrix that acts on the matrix, the solution vector and the
>> righthand side of the problem.
>>
>> The core of the problem I have is that the pressure or potential
>> difference should not be between the leftmost and rightmost site in the
>> physical domain (a standard Dirichlet BC), but between the left- or
>> rightmost site in the physical domain and the corresponding halo site at
>> the opposite side of the system. It should be possible to do this if the
>> entries of the transformation matrix that act on the halo sites can be
>> accessed and modified.
>>
>> Is anything like this possible in PETSc?
>>
>
> The way that periodicity works in PETSc right now for DMDA is that values
> are copied into the halo
> region from the other part of the mesh. Thus, you can just choose to add
> your delta at the right boundary
> and not at the left. The mechanics would be the same as now. Does that
> make sense?
>
>   Thanks,
>
>      Matt
>
>
>> Best regards and many thanks,
>> Oliver
>>
>> --
>> Dr Oliver Henrich
>> Edinburgh Parallel Computing Centre
>> School of Physics and Astronomy
>> University of Edinburgh
>> King's Buildings, JCMB
>> Edinburgh EH9 3FD
>> United Kingdom
>>
>> Tel: +44 (0)131 650 5818
>> Fax: +44 (0)131 650 6555
>>
>> --
>> The University of Edinburgh is a charitable body, registered in
>> Scotland, with registration number SC005336
>>
>>
>>
>>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which their
> -- Norbert Wiener
>
>
> --
> Dr Oliver Henrich
> Edinburgh Parallel Computing Centre
> School of Physics and Astronomy
> University of Edinburgh
> King's Buildings, JCMB
> Edinburgh EH9 3FD
> United Kingdom
>
> Tel: +44 (0)131 650 5818
> Fax: +44 (0)131 650 6555
>
> --
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336
>
>
>
>

--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their