Periodic boundary conditions plus constant

I need to implement boundary conditions in 2D that are periodic plus a constant, i.e.:

u(y + e_1) = u(y) + const1
u(y + e_2) = u(y) + const2

where e_1 and e_2 are the standard unit vectors (1,0) and (0,1), u is a function, and y is in a 2d unit square. However, I know periodic boundary conditions are usually implemented using a map function such as

def map(self, x, y):
        y[0] = x[0] - 1.0
        y[1] = x[1]

but here, y and x are coordinates rather than function values. I haven’t found a mechanism or example that would allow me to do this. Any suggestions or ideas would be appreciated.

Could you split the solution into two parts, e.g.

u(x, y) = u_0(x, y) + u_1(x, y)

where u_0(x, y) satisfies the periodic BC

u_0(0, y) = u_1(1, y)

and u_1(x, y) satisfies the constant offset condition

u_1(0, y) = 0, \\ u_1(1, y) = C

This would give

\begin{align} u(0, y) =& u_0(0, y) \\ u(1, y) =& u_0(1, y) + C = u_0(0, y) + C = u(0, y) + C \end{align}
2 Likes

Interesting idea, but the conditions are not actually separable:
the two constants are not the same, so the conditions of

u(right boundary) = c_1 and u(top boundary) = c_2 are not possible in a continuous manner. The periodicity allows for variation along the right and top boundaries, i.e. more flexibility than the Dirichlet condition allows.

I think I need to access the function values withing the periodic boundary condition, but I’m not aware of a way to access function values.

Gotcha. Are you familiar with @dokken’s dolfinx_mpc? I don’t think this type of boundary condition is currently possible in that library, but it might only require a modest amount of effort to add that capability, e.g. in https://github.com/jorgensd/dolfinx_mpc/blob/main/python/dolfinx_mpc/multipointconstraint.py#L160. You could make a feature request.

Just out of curiosity, would you mind sharing more about the problem you’re solving? I’ve never encountered this type of bc in my own work and I’m curious about how/where it arises. :slight_smile:

Thanks for the suggestion- I haven’t heard of that extension. I’ll make a feature request and in the meantime use a hacky solution like interpolating function values on the other boundary at each timestep and manually feeding them plus a constant in to the periodic boundary as Dirichlet.

Sure- I’m working on a multiscale time-dependent constitutive PDE. In its full multiscale glory, the PDE is
-\nabla \cdot(E^{\epsilon} \nabla u^{\epsilon} + \nu^{\epsilon} \partial t \nabla u^{\epsilon}) = f, x \in \Omega
u^{\epsilon} = 0 , x \in \partial \Omega
u^\epsilon(0) = u^{*} (Initial condition)
where y = \frac{x}{\epsilon} is the “fine scale” variable for \epsilon \ll 1. The idea is that the PDE parameters E^{\epsilon} and \nu^{\epsilon} have variation on scales that are too small to resolve on the full domain \Omega, but they are periodic in y. What we are actually interested in is a map from average strain \int_{\mathbb{T}^d} \nabla u^{\epsilon} \; dy to average stress \int_{\mathbb{T}^d} (E^{\epsilon} + \nu^{\epsilon}\partial t)\nabla u^{\epsilon} \; dy, where both of these are trajectories in time that no longer depend on the small scale variable. This is related to some key terms like “homogenization” and “effective properties.” If you’re unfamiliar with this, I recommend Chapter 12 of “Multiscale Methods” by Pavliotis and Stuart.

Long story short, the way to solve for the average stress given an average strain trajectory involves a PDE at each timestep that takes the following form (the problem is now posed on the microscale, so \epsilon notation is dropped and the forcing, which is assumed to be on a larger scale, vanishes)
-\nabla \cdot ((E + \nu \partial t)\nabla u) = 0 such that u(y + e_i) = u(y) + \overline{\epsilon}_i where \overline{\epsilon} is the average strain at that time. From this one can see that an attempt to separate the problem into one with a periodic boundary condition and one with a Dirichlet boundary condition would just result in a trivial zero solution for the periodic case and constant edges in the Dirichlet case.

I haven’t been able to explain in real detail here, but if you’re really interested, we have some related papers: links may be found here and here.

1 Like