Skip to main content

Section 4.12 Determining Internal Node Values

As seen in the previous section, a finite element solution of a boundary value problem boils down to finding the best values of the constants \(\{C_{j}\}_{j=1}^{n}\text{,}\) which are the values of the solution at the nodes. The interior nodes values are determined by variational principles. Variational principles usually amount to minimizing internal energy. It is a physical principle that systems seek to be in a state of minimal energy and this principle is used to find the internal node values.

Subsection 4.12.1 Variational Principles

For the differential equations that describe many physical systems, the internal energy of the system is an integral. For instance, for the steady state heat equation
\begin{equation} u_{xx}(x,y) + u_{yy}(x,y) = g(x,y)\tag{4.12.1} \end{equation}
the internal energy is the integral
\begin{equation} I[u] = \iint_{R} u_{x}^{2}(x,y) + u_{y}^{2}(x,y) + 2g(x,y) u(x,y) \,dA\text{,}\tag{4.12.2} \end{equation}
where \(R\) is the region on which we are working. It can be shown that \(u(x,u)\) is a solution of (4.12.1) if and only if it is minimizer of \(I[u]\) in (4.12.2).

Subsection 4.12.2 The finite element solution

Recall that a finite element solution is a linear combination of finite element functions:
\begin{equation*} U(x,y) = \sum_{j=1}^{n} C_{j} \Phi_{j}(x,y)\text{,} \end{equation*}
where \(n\) is the number of nodes. To obtain the values at the internal nodes, we will plug \(U(x,y)\) into the energy integral and minimize. That is, we find the minimum of
\begin{equation*} I[U] \end{equation*}
for all choices of \(\{C_{j}\}_{j=1}^{m}\text{,}\) where \(m\) is the number of internal nodes. In this as with any other minimization problem, the way to find a possible minimum is to differentiate the quantity with respect to the variables and set the results to zero. In this case the free variables are \(\{C_{j}\}_{j=1}^{m}\text{.}\) Thus to find the minimizer we should try to solve
\begin{equation*} \frac{\partial I[U]}{\partial C_{j}}= 0 \qquad \textrm{for}\quad 1 \le j \le m\text{.} \end{equation*}
We call this set of equations the internal node equations. At this point we should ask whether the internal node equations can be solved, and if so, is the solution actually a minimizer (and not a maximizer). The following two facts answer these questions. These facts make the finite element method practical:
  • For most applications the internal node equations are linear.
  • For most applications the internal node equations give a minimizer.
We can demonstrate the first fact using an example.

Subsection 4.12.3 Application to the steady state heat equation

If we plug the candidate finite element solution \(U(x,y)\) into the energy integral for the heat equation (4.12.2), we obtain
\begin{equation*} I[U] = \iint_{R} U_{x}(x,y)^{2} + U_{y}(x,y)^{2} + 2g(x,y) U(x,y) \,dA\text{.} \end{equation*}
Differentiating with respect to \(C_{j}\) we obtain the internal node equations
\begin{equation} 0 = \iint_{R} 2 U_{x} \frac{\partial U_{x}}{\partial C_{j}}+ 2 U_{y} \frac{\partial U_{y}}{\partial C_{j}}+ 2 g(x,y) \frac{\partial U}{\partial C_{j}}\,dA \quad\text{for}\quad 1 \le j \le m\text{.}\tag{4.12.3} \end{equation}
Now we have several simplifications. First note that since
\begin{equation*} U(x,y) = \sum_{j=1}^{n} C_{j} \Phi_{j}(x,y)\text{,} \end{equation*}
we have
\begin{equation*} \frac{\partial U}{\partial C_{j}}= \Phi_{j}(x,y)\text{.} \end{equation*}
Also note that
\begin{equation*} U_{x}(x,y) = \sum_{j=1}^{n} C_{j} \frac{\partial}{\partial x}\Phi_{j}(x,y)\text{,} \end{equation*}
and so
\begin{equation*} \frac{\partial U_{x}}{\partial C_{j}}= (\Phi_{j})_{x}\text{.} \end{equation*}
Similarly, \(\frac{\partial U_{y}}{\partial C_{j}}= (\Phi_{j})_{y}\text{.}\) The integral (4.12.3) then becomes
\begin{equation*} 0 = 2 \iint U_{x} (\Phi_{j})_{x} + U_{y} (\Phi_{j})_{y} + g(x,y) \Phi_{j}(x,y) \,dA \quad\text{for}\quad 1 \le j \le m\text{.} \end{equation*}
Next we use the fact that the region \(R\) is subdivided into triangles \(\{T_{i}\}_{i=1}^{p}\) and the functions in question have different definitions on each triangle. The integral then is a sum of the integrals:
\begin{equation*} 0 = 2 \sum_{i=1}^{p} \iint_{T_i}U_{x} (\Phi_{j})_{x} + U_{y} (\Phi_{j})_{y} + g(x,y) \Phi_{j}(x,y) \,dA \quad\text{for}\quad 1 \le j \le m\text{.} \end{equation*}
Now note that the function \(\Phi_{j}(x,y)\) is linear on triangle \(T_{i}\) and so
\begin{equation*} \Phi_{ij}(x,y) = \Phi_{j}|_{T_i}(x,y) = a_{ij}+ b_{ij}x + c_{ij}y\text{.} \end{equation*}
This gives us the simplifications
\begin{equation*} (\Phi_{ij})_{x}(x,y) = b_{ij}\quad \textrm{ and }\quad (\Phi_{ij})_{y}(x,y) = c_{ij}\text{.} \end{equation*}
Also, \(U_{x}\) and \(U_{y}\) restricted to \(T_{i}\) have the form
\begin{equation*} U_{x} = \sum_{k=1}^{n} C_{k} b_{ik}\quad \textrm{ and }\quad U_{y} = \sum_{k=1}^{n} C_{k} c_{ik}\text{.} \end{equation*}
The internal node equations then reduce to
\begin{equation*} 0 = \sum_{i=1}^{p} \iint_{T_i}\left( \sum_{k=1}^{n} C_{k} b_{ik}\right) b_{ij}+ \left( \sum_{k=1}^{n} C_{k} c_{ik}\right) c_{ij}+ g(x,y) \Phi_{ij}(x,y) \,dA \quad\text{for}\quad 1 \le j \le m\text{.} \end{equation*}
Now notice that \(\left( \sum_{k=1}^{n} C_{k} b_{ik}\right) b_{ij}\) is just a constant on \(T_{i}\text{,}\) and, thus, we have
\begin{equation*} \iint_{T_i}\left( \sum_{k=1}^{n} C_{k} b_{ik}\right) b_{ij}+ \left( \sum_{k=1}^{n} C_{k} c_{ik}\right) c_{ij}= \left[ \left( \sum_{k=1}^{n} C_{k} b_{ik}\right) b_{ij}+ \left( \sum_{k=1}^{n} C_{k} c_{ik}\right) c_{ij}\right] A_{i}\text{,} \end{equation*}
where \(A_{i}\) is just the area of \(T_{i}\text{.}\) Finally, we apply the Three Corners rule to make an approximation to the integral
\begin{equation*} \iint_{T_i}g(x,y) \Phi_{ij}(x,y) \,dA\text{.} \end{equation*}
Since \(\Phi_{ij}(x_{k},y_{k})=0\) if \(k\not=j\) and even \(\Phi_{ij}(x_{j},y_{j})=0\) if \(T_{i}\) does not have a corner at \((x_{j},y_{j})\text{,}\) we get the approximation
\begin{equation*} \Phi_{ij}(x_{j},y_{j})g(x_{j},y_{j}){A_i}/{3}\text{.} \end{equation*}
If \(T_{i}\) does have a corner at \((x_{j},y_{j})\) then \(\Phi_{ij}(x_{j},y_{j})=1\text{.}\)
Summarizing, the internal node equations are
\begin{equation*} 0 = \sum_{i=1}^{p} \left[ \left( \sum_{k=1}^{n} C_{k} b_{ik}\right) b_{ij}+ \left( \sum_{k=1}^{n} C_{k} c_{ik}\right) c_{ij}+ \frac{1}{3}g(x_{j},y_{j})\Phi_{ij}(x_{j},y_{j}) \right] A_{i} \quad\text{for}\quad 1 \le j \le m\text{.} \end{equation*}
While not pretty, these equations are in fact linear in the unknowns \(\{C_{j}\}\text{.}\)

Subsection 4.12.4 Experiment

Download the program myfiniteelem.m (Program A.2.3). This program produces a finite element solution for the steady state heat equation without source term:
\begin{equation*} u_{xx}+ u_{yy}= 0\text{.} \end{equation*}
To use it, you first need to set up the region and boundary values by running a script such as mywasher.m (Program A.2.8) or mywedge.m (Program A.2.9). Try different settings for the boundary values z. You will see that the program works no matter what you choose.

Exercises 4.12.5 Exercises

For this exercise, use the region that you made in Exercise 4.11.6.1.
Does the second plot look like it satisfies a variational principle? Why or why not? Turn in the two plots and your answers to the questions.