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?


2 Answers 2


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}$


$\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


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


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.

  • $\begingroup$ 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}$. $\endgroup$
    – andybauer
    Jan 22, 2013 at 18:55
  • $\begingroup$ 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. $\endgroup$ Jan 23, 2013 at 5:39
  • $\begingroup$ 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. $\endgroup$
    – andybauer
    Jan 23, 2013 at 15:44

To add up to Nathan's great answer with the variational reasoning, one often needs algorithmic details when implementing finite elements. For example,

Algorithm 1

I also have a more detailed explanation on the subject in my personal notes. Please see the chapter "Constrained Linear Systems".


Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.