To find the finite element solution, we can either start with the
strong form and derive the weak form, or we can start with a weak form
derived from a variational principle.

Let us assume that the approximate solution is $\mathbf {u} _{h}(\mathbf {x} )$ and plug
it into the ODE. We get

$AE{\cfrac {d^{2}\mathbf {u} _{h}}{dx^{2}}}+a\mathbf {x} =R_{h}(\mathbf {x} )$ where $R_{h}$ is the residual . We now try to minimize the residual in a weighted average sense

$\int _{0}^{L}R_{h}(\mathbf {x} )\mathbf {w} (\mathbf {x} )~dx=0$ where $\mathbf {w} (\mathbf {x} )$ is a weighting function. Notice that this equation is similar to equation (5) (see 'Weak form: integral equation') with $\mathbf {w}$ in place of the variation $\delta \mathbf {u}$ . For the two equations to be equivalent, the weighting function must also be such that $\mathbf {w} (0)=0$ .

Therefore the approximate weak form can be written as

${\int _{0}^{L}AE{\cfrac {d\mathbf {u} _{h}}{dx}}{\cfrac {d\mathbf {w} }{dx}}~dx=\int _{0}^{L}\mathbf {q} \mathbf {w} ~dx+\left.{\boldsymbol {R}}~\mathbf {w} \right|_{x=L}~.}$ In Galerkin's method we assume that the approximate solution can
be expressed as

$\mathbf {u} _{h}(\mathbf {x} )=a_{1}\varphi _{1}(\mathbf {x} )+a_{2}\varphi _{2}(\mathbf {x} )+\dots +a_{n}\varphi _{n}(\mathbf {x} )=\sum _{i=1}^{n}a_{i}\varphi _{i}(\mathbf {x} )~.$ In the Bubnov-Galerkin method , the weighting function is chosen to be of the same form as the approximate solution (but with arbitrary coefficients),

$\mathbf {w} (\mathbf {x} )=b_{1}\varphi _{1}(\mathbf {x} )+b_{2}\varphi _{2}(\mathbf {x} )+\dots +b_{n}\varphi _{n}(\mathbf {x} )=\sum _{j=1}^{n}b_{j}\varphi _{j}(\mathbf {x} )~.$ If we plug the approximate solution and the weighting functions into
the approximate weak form, we get

$\int _{0}^{L}AE\left(\sum _{i=1}^{n}a_{i}{\cfrac {d\varphi _{i}}{dx}}\right)\left(\sum _{j=1}^{n}b_{j}{\cfrac {d\varphi _{j}}{dx}}\right)~dx=\int _{0}^{L}\mathbf {q} \left(\sum _{j=1}^{n}b_{j}\varphi _{j}\right)~dx+\left.{\boldsymbol {R}}~\left(\sum _{j=1}^{n}b_{j}\varphi _{j}\right)\right|_{x=L}~.$ This equation can be rewritten as

$\sum _{j=1}^{n}b_{j}\left[\int _{0}^{L}AE\left(\sum _{i=1}^{n}a_{i}{\cfrac {d\varphi _{i}}{dx}}{\cfrac {d\varphi _{j}}{dx}}\right)~dx\right]=\sum _{j=1}^{n}b_{j}\left[\int _{0}^{L}\mathbf {q} \varphi _{j}~dx+\left.\left({\boldsymbol {R}}~\varphi _{j}\right)\right|_{x=L}\right]~.$ From the above, since $b_{j}$ is arbitrary, we have

$\int _{0}^{L}AE\left(\sum _{i=1}^{n}a_{i}{\cfrac {d\varphi _{i}}{dx}}{\cfrac {d\varphi _{j}}{dx}}\right)~dx=\int _{0}^{L}\mathbf {q} \varphi _{j}~dx+\left.{\boldsymbol {R}}~\varphi _{j}\right|_{x=L}~,~j=1\dots n.$ After reorganizing, we get

$\sum _{i=1}^{n}\left[\int _{0}^{L}{\cfrac {d\varphi _{j}}{dx}}AE{\cfrac {d\varphi _{i}}{dx}}~dx\right]a_{i}=\int _{0}^{L}\varphi _{j}\mathbf {q} ~dx+\left.\varphi _{j}{\boldsymbol {R}}\right|_{x=L}~,~j=1\dots n$ which is a system of $n$ equations that can be solved for the unknown coefficients $a_{i}$ . Once we know the $a_{i}$ s, we can use them to compute approximate solution. The above equation can be written in matrix form as

$\mathbf {K} \mathbf {a} =\mathbf {f} \qquad \leftrightarrow \qquad K_{ji}a_{i}=f_{j}$ where

$\mathbf {K} =\int _{0}^{L}\mathbf {B} ^{T}\mathbf {D} \mathbf {B} ~dx\qquad \leftrightarrow \qquad K_{ji}=\int _{0}^{L}{\cfrac {d\varphi _{j}}{dx}}AE{\cfrac {d\varphi _{i}}{dx}}~dx$ and

$f_{j}=\int _{0}^{L}\varphi _{j}\mathbf {q} ~dx+\left.\varphi _{j}{\boldsymbol {R}}\right|_{x=L}~.$ The problem with the general form of the Galerkin method is that the
functions $\varphi _{i}$ are difficult to determine for complex domains.