# How to get the linear equation system for finite element method from the variational formulation

Let the problem be $$-u'' + a(x) u = f , \;x \in \Omega = ]0,1[ , \;u(0) = \alpha ; \;u(1) = \beta,$$ where $f \in L^2(\Omega) , a(x) \geq a_0 > 0 , a(x) \in L^{\infty}(\Omega).$

This problem admits a unique solution in $V = H^1(\Omega).$

My question is:

How we can use the finite elements $\mathbb{P}_1$ to prove that this problem, in the finite element space, is equivalent to a linear system $A U = b$ .

BTW: for this, we can use the relation for linear function $\psi$: $$\dfrac{1}{h} \int_{x_i}^{x_{i+1}} \psi(x) dx = \psi\left(\dfrac{x_i + x_{i+1}}{2}\right).$$

-
You just write the 2nd derivative as a second-order scheme: $$u''(x)\approx \frac{u(x+\Delta) + u(x-\Delta) - 2u(x)}{2\Delta}$$ for the uniform grid on $[0,1]$. The only thing you would need to determine this system of linear equations are boundary conditions at $x = 0$ and $x = 1$ which you do have. What's the problem? –  Ilya Jun 1 '13 at 12:34
i want a finite element NOT difference finies!! –  jijiii Jun 1 '13 at 16:53
Here how we can use finite element? –  jijiii Jun 1 '13 at 16:53
Do you know what is the nodal basis function for $\mathbb{P}^1$ elements? –  Shuhao Cao Jun 1 '13 at 18:17
yes. $$V_h = \{ v_h \in H^1: v_h \in \mathbb{P}^1, for all x \in I_i\}$$ $I_i = [x_i,x_{i+1}], i=0,...,n.$ and the basic of $V_h$ is the functions $w_j$ such that $w_j(x_i)=1$ if $i=j$ and $w_j(x)=0$ if $i \neq j.$ My problem is to compute the matrix A. –  jijiii Jun 1 '13 at 19:23

The nodal basis function for $\mathbb{P}_1$-elements is exactly as you said in the comments, a.k.a. hat functions (see Figure 11.2 here): $$w_j(x_i) = \begin{cases}1\quad \text{ when }i = j, \\ 0\quad \text{ otherwise },\end{cases}$$ and the support of a hat function $w_j$ is $(x_{j-1},x_{j+1})$ in which the hat function is linear. First we assume the finite element solution $$u_h(x) = \alpha w_0(x) + \beta w_n(x) + \sum^{n-1}_{n=1} u_i w_i(x),\tag{1}$$ so that $u_h$ satisfies the boundary condition, and $u_i$'s above are coefficients (degrees of freedom sometimes we call them) we want to find. Now we wanna plugging the hat function one by one in the weak formulation: $$\int_0^1 \left(u'(x)v'(x) + a(x)u(x)v(x)\right)\,dx = \int_0^1 f(x) v(x)\,dx.\tag{2}$$ Plugging (1) into (2), then replacing $v(x)$ by hat functions $w_j$, $j = 1,\ldots,n-1$ (do not include the hat functions on the end points for the test function space has zero boundary condition), we have: $$\int_0^1 \left(\big(\sum^{n-1}_{n=1} u_i w_i(x)\big)'w_j'(x) + a(x)\big(\sum^{n-1}_{n=1} u_i w_i(x)\big)w_j(x)\right)\,dx = \int_0^1 f(x) w_j(x)\,dx,$$ pulling the summation outside the integral for $u_i$'s are numbers: $$\sum^{n-1}_{n=1} u_i \int_0^1 \left(w'_i(x)w_j'(x) + a(x) w_i(x) w_j(x)\right)\,dx = \int_0^1 f(x) w_j(x)\,dx,$$ Now we have a linear system: \begin{gather} &AU = b, \quad \text{ where } \\ &A_{ji} = \int_0^1 \left(w'_i(x)w_j'(x) + a(x) w_i(x) w_j(x)\right)\,dx : \text{ an } (n-1)\times (n-1) \text{ matrix } , \\ &U_{i} = u_i: \text{ an } (n-1)\text{-column vector } , \\ &b_{j} = \int_0^1 f(x) w_j(x)\,dx: \text{ an } (n-1)\text{-column vector } . \end{gather} Now the last task is to evaluate: $$A_{ji} = \int_0^1 \left(w'_i(x)w_j'(x) + a(x) w_i(x) w_j(x)\right)\,dx,$$ and $A_{ji}$ is zero if two hat functions do not overlap. $w_i(x)$ and $w_j(x)$ overlap when $i= j,j-1,j+1$, so the rest work is evaluate: \begin{aligned} &A_{ii} = \int_{x_{i-1}}^{x_{i+1}} \left(w'_i(x)w_i'(x) + a(x) w_i(x) w_i(x)\right)\,dx, \\ &A_{i,i-1} = \int_{x_{i-2}}^{x_{i+1}} \left(w'_{i-1}(x)w_i'(x) + a(x) w_{i-1}(x) w_i(x)\right)\,dx, \\ &A_{i+1,i} = \int_{x_{i-1}}^{x_{i+2}} \left(w'_{i}(x)w_{i+1}'(x) + a(x) w_{i}(x) w_{i+1}(x)\right)\,dx. \end{aligned}\tag{3} I suggest you try this yourself with the integration formula you gave for linear functions.

UPDATE 1: Nodal basis "hat" functions are: $$w_i(x)=\frac{x−x_{i−1}}{h} \;\text{ in } [x_{i−1},x_i], \quad \text{and}\quad w_i(x)=−\frac{x−x_{i+1}}{h} \;\text{ in } [x_i,x_{i+1}].$$ First integration in $A_{ii}$, notice here we assume all $x_i$'s are evenly spaced, $x_i - x_{i-1} = x_{i+1}-x_i = h$: $$\int_{x_{i-1}}^{x_{i+1}} w'_i(x)^2 = \frac{x_i−x_{i−1}}{h^2} + \frac{x_{i+1}−x_{i}}{h^2} = \frac{2}{h}.$$ For the second integration, we can't do much, unless we assume $a$ is a piecewise constant, otherwise we have to use quadrature rule, and you can't pull $w_i^2$ out. Assume $a= a_{i-1}$ in $[x_{i−1},x_i]$, and $a = i+1$ in $[x_i,x_{i+1}]$ $$\int_{x_{i-1}}^{x_{i+1}}a(x) w_i(x)^2 = \frac{a_{i-1}h}{3} + \frac{a_i h}{3} = \frac{h(a_{i-1}+a_i)}{3}.$$ For $A_{i,i-1}$ and $A_{i,i+1}$, beware of the overlapping part of different hat functions, so $A_{i,i-1}$ is really the integration from $x_{i-1}$ to $x_i$, and $A_{i,i+1}$ is really the integration from $x_{i}$ to $x_{i+1}$ .

UPDATE 2: The numerical integration of an arbitrary polynomial on interval $[\alpha,\beta]$ can use the following handy rules for a single term: $$\int^{\alpha}_{\beta} w_{\alpha}(x)^m w_{\beta}(x)^n\,dx = \frac{m!\,n!}{(m+n+1)!} |\beta-\alpha|,$$ where $w_a(x)$ is 1 at $a$ and 0 at $b$ and linear within, similarly $w_b$ is the hat function for $b$, in your case, if coefficient is constant: $$\int_{x_{i-1}}^{x_{i}} w_i(x)^2\,dx = \int_{x_{i-1}}^{x_{i}} w_i(x)^2 w_{i-1}(x)^0\,dx= \frac{2!\,0!}{(2+0+1)!} |x_{i} -x_{i-1} | = \frac{1}{3}h,$$ and $$\int_{x_{i-1}}^{x_{i}} w_i(x) w_{i-1}(x) \,dx= \frac{1!\,1!}{(1+1+1)!} |x_{i} -x_{i-1} | = \frac{1}{6}h.$$

-
Okay, we have $$A_{ii} = \int_{x_{i-1}}^{x_i} [w'_i(x)]^2 dx + \int_{x_i}^{x_{i+1}} [w'_i(x)]2 dx + \int_{x_{i-1}}^{x_{i}} a(x) [w_i(x)]^2 dx + \int_{x_i}^{x_{i+1}} a(x) [w_i(x)]^2 dx$$ $$= \int_{x_{i-1}}^{x_i} \dfrac{1}{h^2} dx + \int_{x_i}^{x_{i+1}} \dfrac{1}{h^2} dx + \int_{x_{i-1}}^{x_i} a(x) \dfrac{1}{h^2} dx + \int_{x_i}^{x_{i+1}} a(x) \dfrac{1}{h^2} dx$$ $$= \dfrac{2}{h} + \dfrac{1}{h^2} \int_{x_{i-1}}^{x_{i+1}} a(x) dx.$$ But who is $A_ii?$ what we can do with $a(x)$ and how we can use the identité $$\dfrac{1}{h} \int_{x_i}^{x_{i+1}} \psi (x) dx = \psi(\dfrac{x_i + x_{i+1}}{2})?$$ –  jijiii Jun 1 '13 at 22:16
@jijiii The integration formula you gave only works for linear function or constant function, you can leave $\displaystyle \int_{x_{i-1}}^{x_{i+1}} a(x) w_i^2(x)$ as it were. BTW did you read my comments about accepting answers? –  Shuhao Cao Jun 1 '13 at 22:19
yes, and i accept the answer. But i don't understand your last answer. –  jijiii Jun 1 '13 at 22:24
hit: i consider $w_i(x) = \dfrac{x-x_i}{h}$ in $[x_{i-1},x_i]$ and $w_i(x) = - \dfrac{x-x_i}{h}$ in $[x_i,x_{i+1}]$ and 0 if $x$ is not in $[x_{i-1},x_i]$ or $[x_{i} , x_{i+1}]$ –  jijiii Jun 1 '13 at 22:30
@jijiii Should be $w_{i}(x) = \dfrac{x-x_{i-1}}{h}$ in $[x_{i-1},x_{i}]$, and $w_{i}(x) = -\dfrac{x -x_{i+1}}{h}$ in $[x_{i},x_{i+1}]$, for $w_i(x) = 1$ at $x_i$. –  Shuhao Cao Jun 1 '13 at 22:47

i have a question for the solution of this exercice. How we found that $$A_{i,i-1}=\int_{x_{i-2}}^{x_{i+1}} (w'_{i-1}(x) w'_i(x) + a(x) w_{i-1}(x) w_i(x) dx$$? and why we should found this expression of $A_{i,i-1}$ thank's for the help.

-