0
$\begingroup$

Using FEM, I want to approximate the Laplacian

$$u = \nabla \cdot \nabla h \, ,$$

where $h(x,y)$ is an FEM approximated scalar field on the same mesh, i.e. piecewise differentiable.

I am using MOOSE to solve the following (hopefully correct!) weak formulation of the above equation:

$$\int_\Omega uw_i \, \text{d}\Omega = -\int_\Omega \nabla h \cdot \nabla w_i \, \text{d}\Omega + \int_\Gamma \nabla h \cdot \vec{n}w_i \, \text{d}\Gamma \,,$$

where $w_i$ is the i-th test function and $\vec{n}$ is the normal vector on the boundary $\Gamma$ of the computational domain $\Omega$. I am not using any boundary conditions on $u$.

The FEM solution that I obtain for this weak formulation is inaccurate at the boundaries (see below). I have tested this by prescribing the field $h=\text{sin}(x)$ on a 2D rectangular domain from $-2\pi \le x \le 2\pi$ and solving with quadratic Lagrange shape functions for both $h$ and $u$. (it is also inaccurate using linear shape functions, but it is easier to spot with quadratic shape functions).

Is it possible to obtain an accurate $u$ at the boundaries with this approach? Am I missing a boundary condition/contribution?

h(x) = sin(x) FEM solution for u


Background:

I am dealing with a thin film approximation of the Navier-Stokes equations in 2D FEM, where the film height is modeled by a dependent variable $h(x,y)$. In order to include surface tension as a contribution to the pressure gradient, I want to find the curvature of the film surface $u(x,y)$ that I approximate by the second derivative

$u \approx \nabla \cdot \nabla h$ .

My idea is to introduce curvature $u$ as a variable, solve the above equation by FEM, and couple $u$ into the momentum equation.

$\endgroup$
9
  • $\begingroup$ Let me see if I understand what you want to achieve. You have a function $h(x, y)$ and you want to compute it's Laplacian over a domain? $\endgroup$ – nicoguaro Mar 29 at 18:04
  • $\begingroup$ Are you ignoring the boundary term from your weak form? The formulation you have should show up with a Neumann boundary term for GFEM. $\endgroup$ – wwfe Mar 29 at 18:11
  • $\begingroup$ @nicoguaro Yes, but let me reformulate: I want to calculate u by FEM so that it approximates the Laplacian of h(x, y) ("source term" on the right hand side of the equation). I don't have a closed form of h(x, y) though, because it is actually a variable in another FEM problem. $\endgroup$ – pimpom Mar 29 at 23:10
  • $\begingroup$ @wwfe Can you tell me what boundary term you mean and what it should look like? I don't have a boundary term with u since I am not integrating by parts on the LHS of the pde. $\endgroup$ – pimpom Mar 29 at 23:13
  • $\begingroup$ @pimpom the boundary represented by the integral over $\Gamma$ is in your equation, just because u isn't there doesn't mean it isn't a boundary term, it is still part of the residual. $\endgroup$ – wwfe Mar 30 at 3:04

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.