Implementing both Dirichlet and Neumann boundary conditions for a fourth-order PDE

Hello,
Consider a fourth-order PDE in a unit square \Omega described by cartesian coordinates x^1, x^2 for the function f
\nabla^2 H\left(f, {\bf \nabla } f, {\bf \nabla }^2 f\right) + \cdots = 0,
where H is a non-linear function of f and of its derivatives up to second order.

Given that the PDE is of fourth order, I have two boundary conditions (BCs) on \partial \Omega:

  • a Dirichlet BC f({\bf x}) = h({\bf x}) on \partial \Omega
  • a Neumann BC n^i \partial_i H = k({\bf x}) on \partial \Omega, where \bf n is the unit normal to \partial \Omega.

In the variational formulation, I multiply the PDE by a test function \nu, integrate, and obtain a functional

F = \langle (\nabla^2 H) q \rangle + \cdots. I integrate by parts, and get
F = - \langle (\partial_i H) \partial_i q \rangle_{\Omega} + \langle n ^i (\partial_i H) q \rangle_{\partial \Omega} + \cdots.

What is the correct way to impose the Neumann BCs here ? Imposing them by setting \langle n ^i (\partial_i H) q \rangle_{\partial \Omega} \rightarrow \langle k\, q \rangle_{\partial \Omega} would have no effect, because the Dirichlet BCs imply that q=0 on \partial \Omega.

Thank you

Two methods spring to mind:

  1. Employ an H^2 conforming FE basis which preserves C^1 continuity of your approximation. Then you can impose both boundary conditions in a strong sense. I’m not sure of the implementation status of these bases from Basix and their support with dolfinx. In the past people have been referred to @kamensky’s tIGAr project, e.g., this demo; however, this project has not been updated for a long time.
  2. Employ weak imposition of C^1 continuity via, e.g., Nitsche’s method. There’s a demo, here.