I am trying to solve simple scalar biharmonic equation using bubnov-galerkin finite element method. I am using $H^2$ conforming basis functions. I was wondering that if anyone can give me some pointers on how can I further debug my algorithm.I proceeded by converting the strong form of equation to weak form:\begin{align}\Delta^2 u &= f \\u &= g \;\;\; \text{on} \;\;\; \Gamma_D \\\Gamma_N &= \emptyset\end{align}As an example problem, I am trying to solve:\begin{align}u &= \cos(4 \pi x) \cos(4 \pi y)\end{align}Boundary conditions are implemented using penalty method. The weak form is follows:\begin{align}a(u,v) &= L(v) \\a(u,v) &= \sum_{\Omega_K \in \mathcal T} \int_{\Omega_K} \Delta \psi_{i} \Delta \psi_{j} + \sum_{E \in \mathcal E} \gamma \int_{E} \psi_{i} \psi_{h,j} \\L(v) &= \sum_{\Omega_K \in \mathcal T} \int_{\Omega_K} \psi_{i} f + \sum_{E \in \mathcal E} \gamma \int_{E} \psi_{i} g\\\end{align}Problem: When I solve this equation on a square domain, my L2 error barely converges (convergence rate ~ 0.3). However, if I solve $u = \sin(4 \pi x) \sin(4 \pi y)$, I get correct convergence rates. I have tried the following things to debug my code:I solves the poisson equation, so replaced the first integrand with the stiffness integrand. I get correct convergence rates. I concluded from this that my penalty method implementation is rightI am using Cartesian grid, so Jacobian lines up with expectationsSince, my shape functions are defined in parametric ($\xi,\eta$) space, I had to work out the transformation for the laplacian. I tested this transformation on polar coordinates.Many thanks.Apologies. Editted.
Boundary conditions in conforming Galerkin method for biharmonic equation
pde;finite element;boundary conditions