Periodicity and Robin bcs

General
1

Hi, I am trying to setup a periodic bc on a simple geometry, a 2d square, where 2 opposite edges have a periodic condition. In one of the other 2 edges of the square, I apply a Robin bc, thus 2 corner nodes are subjected to both periodic condition and Robin bc. If I try to solve the Laplace problem, I see that the solution is wrong. If instead I remove the corner node from the periodicity indicator lambda, the solution is what I expect, i.e. an equipotential solution = 1 in the whole domain.

from typing import Any
import dolfinx
import dolfinx_mpc
import ufl
from dolfinx.fem import functionspace
from dolfinx_mpc import MultiPointConstraint
from mpi4py import MPI
from utils import FenicsxPlotter

mesh1 = dolfinx.mesh.create_unit_square(MPI.COMM_WORLD, 10, 10)
def periodic_slave(x: Any) -> Any:
    return np.logical_and(np.isclose(x[1], 1),True)        #this gives wrong solution
    #return np.logical_and(np.isclose(x[1], 1),x[0]>1e-5)    #this is ok
def mapping(slave_x: Any) -> Any:
    master_x = np.copy(slave_x)
    master_x[0] = 0.0
    return master_x
V = functionspace(mesh1, ("Lagrange", 1))
mpc = MultiPointConstraint(V)
u = ufl.TrialFunction(V)
v = ufl.TestFunction(V)
facets = dolfinx.mesh.locate_entities(mesh1, 1, lambda x: np.isclose(x[1], 0.0))
facet_tag = dolfinx.mesh.meshtags(
    mesh1,
    1,
    np.unique(facets),
    1,
)
ds = ufl.Measure("ds", domain=mesh1, subdomain_data=facet_tag)
a = ufl.dot(ufl.grad(u), ufl.grad(v)) * ufl.dx + ufl.inner(u, v) * ds
L = dolfinx.fem.Constant(mesh1, 0.0) * v * ufl.dx + v * ds
mpc.create_periodic_constraint_geometrical(V, periodic_slave, mapping, [])
mpc.finalize()
problem = dolfinx_mpc.LinearProblem(
    a,
    L,
    mpc,
    bcs=[],
    petsc_options={"ksp_type": "preonly", "pc_type": "lu"},
)
uh = problem.solve()

The solution I get is:

wrong
wrong734×632 135 KB

the solution I get excluding the corner node in periodic_slave:
ok
ok891×632 151 KB

I was expecting that only Dirichlet bcs should be excluded when applying periodic conditions, as the signature of create_periodic_constraint_geometrical also suggest. What am I missing?

2

You’re trying to be periodic in the x-direction, correct? Then I think you just made a small typo:

return np.logical_and(np.isclose(x[1], 1),True)

should be

return np.logical_and(np.isclose(x[0], 1),True)
3

The periodicity I wanted to try is along y axis, so I was setting the left boundary (x=0) as master and the right boundary (x=1) as slave boundary. The Robin bc is applyed to the bottom boundary (y=0). Thus

return np.logical_and(np.isclose(x[1], 1),True)

should be ok.

4

x=0 as master and x=1 as slave means periodic in x…?

In any case, that line is testing y=1, which you don’t want

1 Like
5

Ah, sorry, you’re right. I guess I’ve not taken enough coffee this morning :slight_smile: Setting it up correctly, it works as expected.

1 Like