# How to properly apply non-homogeneous Dirichlet boundary conditions with FEM?

In general, Dirichlet boundary conditions won't be satisfied exactly for FEM for non-homogeneous boundary conditions. The FEM codes I've seen set the degrees of freedom to interpolate the Dirichlet boundary condition but I haven't found any mathematical justification for this. It seems to me that setting essential boundary conditions should probably minimize some functional of the error (e.g. minimize $||u -u_h||$ over the portion of the boundary that the Dirichlet BC is applied ) even though this would be more computationally expensive.

Is there any justification for setting the BC like this and if so, what would the proper norm be?

-

There is mathematical justification for setting Dirichlet boundary degrees of freedom to a value. However, you should adjust your variational form accordingly. If you are looking at a general problem, say:

Find $u\in\mathcal{U}$ such that

$a(u,w)=l(w) \ \ \forall w\in\mathcal{V}$

where

$\mathcal{U}=\{u:\int \nabla u^2 < \infty, u=g\text{ on }\Gamma_D\}$

$\mathcal{V}=\{u:\int \nabla u^2 < \infty, u=0\text{ on }\Gamma_D\}$

Instead we can write $u = v + g$ where $v\in\mathcal{V}$ and $g$ is the Dirichlet condition. Then the variational form becomes

$a(v+g,w)=l(w)$

or by using the linearity of $a(.,.)$

$a(v,w)=l(w)-a(g,w)$

In a finite element code, you can form your element stiffness matrix as if there were no boundary conditions. Then you take the column of the local matrix which corresponds to the Dirichlet boundary condition, scale it by the coefficient you want to enforce, and subtract it from the right-hand-side. This is the discrete form of what I wrote above, $-a(g,w)$. Then you zero out that column and the corresponding Dirichlet row, placing a 1 in the diagonal and the coefficient you wish to enforce. This decouples the equation from the system and yet sets the value you wish to enforce.

I recommend The Finite Element Method: Linear Static and Dynamic Finite Element Analysis, by Tom Hughes. He has an expanded discussion of this issue starting on page 8.

-
I'm not sure I follow your manipulations with the local matrix. In the Hughes book the example is for a 1 dimensional problem where it is trivial to satisfy the Dirichlet boundary condition (i.e. equation 1.6.4) but I don't see yet how it should work for 2D or 3D when $g$ can't be satisfied exactly on the boundary by a given implementations shape functions (e.g $g(x) = x^2$ and using piecewise linear shape functions). Also, I know I'm nitpicking by you should probably change to $u=v+g^h$ with $g^h \in \mathcal{U}$. – andybauer Jan 22 '13 at 18:55
Yes, you have the right idea. I only wrote the continuous part above and did not specify discretization. By setting coefficients in the solution, we are assuming that the Dirichlet condition we wish to enforce may be represented by the function space we choose. In your example of $g(x)=x^2$, say on the bottom part of the unit domain, there are two ways I know that this would be handled: (1) Setting the coefficient to whatever $g(x)$ is at the dof location (2) $L^2$-projection of the function on the boundary to get the coefficients. – Nathan Collier Jan 23 '13 at 5:39
Thanks -- I guess what I was trying to get at in my poorly worded question was whether we should do (1) or (2). (1) seems to be the way that I've seen done in the FEM codes I've looked at but (2) seems like it would result in a better approximation. – andybauer Jan 23 '13 at 15:44