Question 2 of 7: Bubnov–Galerkin cubic solution for an axially loaded cantilevered bar

Question 2: Bubnov–Galerkin cubic solution for an axially loaded cantilevered bar (20 marks)

Question text not reproduced: the examination questions are © Engineers and Geoscientists BC. Open the official past paper (linked at the top of this page) to read the question, then follow the worked solution below.

Given. The bar is governed by the second-order equilibrium equation above with one essential and one natural boundary condition, and a cubic trial field is called for.

Find. The cubic Bubnov–Galerkin approximation $u(x)$, together with the resulting tip displacement and axial stress distribution.

Approach. Pick a cubic trial function that satisfies the essential boundary condition identically, use the trial functions themselves as weighting functions (that is what makes the method Bubnov–Galerkin), integrate the weighted residual by parts so the natural boundary condition enters naturally, and solve the resulting $3\times 3$ system.

1. Choose an admissible cubic trial field. Take $$\tilde u(x) = a_1x + a_2x^{2} + a_3x^{3} = \sum_{i=1}^{3}a_i\phi_i,\qquad \phi_i(x) = x^{i}$$ Every term vanishes at $x=0$, so the essential condition $u(0)=0$ is satisfied a priori and no constant term is admitted. The natural condition is left to the weak form; it must not be imposed on the trial function.
2. Write the weighted-residual statement. The residual of the governing equation is $R(x) = EA\,\tilde u''(x) + cx$. In the Bubnov–Galerkin method the weighting functions are the same functions used to build the trial field, $w_i = \phi_i = x^{i}$, and each weighted residual is driven to zero: $$\int_{0}^{L}\phi_i\left[EA\,\tilde u'' + cx\right]dx = 0,\qquad i=1,2,3$$
3. Integrate by parts to admit the natural boundary condition. Applying integration by parts to the second-derivative term, $$\Bigl[EA\,\phi_i\,\tilde u'\Bigr]_{0}^{L} -\int_{0}^{L}EA\,\phi_i'\,\tilde u'\,dx + \int_{0}^{L}cx\,\phi_i\,dx = 0$$ The boundary term vanishes at $x=0$ because $\phi_i(0)=0$, and at $x=L$ it becomes $\phi_i(L)\,EA\,\tilde u'(L) = P\,L^{i}$ by the given natural condition. The weak statement is therefore $$EA\int_{0}^{L}\phi_i'\,\tilde u'\,dx = P\,L^{i} + c\int_{0}^{L}x\,\phi_i\,dx$$
4. Assemble the coefficient system. With $\phi_i' = i\,x^{i-1}$ the integrals are elementary, $K_{ij} = EA\!\int_0^L ij\,x^{i+j-2}dx = EA\,\dfrac{ij\,L^{i+j-1}}{i+j-1}$ and $F_i = PL^{i} + c\,\dfrac{L^{i+2}}{i+2}$, giving $$EA\begin{bmatrix}L & L^{2} & L^{3}\\[2pt] L^{2} & \tfrac{4}{3}L^{3} & \tfrac{3}{2}L^{4}\\[2pt] L^{3} & \tfrac{3}{2}L^{4} & \tfrac{9}{5}L^{5}\end{bmatrix}\begin{Bmatrix}a_1\\a_2\\a_3\end{Bmatrix}=\begin{Bmatrix}PL+\tfrac{1}{3}cL^{3}\\[2pt] PL^{2}+\tfrac{1}{4}cL^{4}\\[2pt] PL^{3}+\tfrac{1}{5}cL^{5}\end{Bmatrix}$$ The matrix is symmetric, as it must be for a self-adjoint operator — that symmetry is the practical signature of the Bubnov choice of weights.
5. Solve for the three coefficients. Solving the system (and cancelling the common $EA$) gives $$a_1 = \frac{P}{EA} + \frac{cL^{2}}{2EA},\qquad a_2 = 0,\qquad a_3 = -\frac{c}{6EA}$$ so that the Bubnov–Galerkin cubic is $$\boxed{\;u(x) = \frac{Px}{EA} + \frac{c}{6EA}\left(3L^{2}x - x^{3}\right)\;}$$ The quadratic coefficient vanishing is a real result, not an oversight: the loading $q=cx$ produces no linear term in the axial force distribution beyond the constant $P$.
6. Check the residual and interpret the accuracy. Substituting back, $R(x) = EA\left(2a_2 + 6a_3x\right) + cx = EA\left(-\tfrac{c}{EA}x\right) + cx = 0$ identically. The residual is not merely orthogonal to the three weighting functions — it is zero everywhere, so the cubic is the exact solution of the boundary value problem. That happens because integrating $EA\,u'' = -cx$ twice produces a cubic, and the trial space already contains it. Galerkin returns the exact answer whenever the exact solution lies in the trial space; otherwise it returns the best approximation in the energy norm.
7. Extract the engineering quantities. Differentiating, the axial strain and stress are $$\varepsilon(x) = u'(x) = \frac{P}{EA} + \frac{c\left(L^{2}-x^{2}\right)}{2EA},\qquad \sigma(x) = E\,\varepsilon(x) = \frac{P}{A} + \frac{c\left(L^{2}-x^{2}\right)}{2A}$$ The stress is largest at the built-in end, $\sigma(0) = P/A + cL^{2}/(2A)$, where the whole distributed load $\int_0^L cx\,dx = cL^{2}/2$ has accumulated, and falls to $\sigma(L) = P/A$ at the free end, which is exactly the imposed natural condition. The tip displacement is $$u(L) = \frac{PL}{EA} + \frac{cL^{3}}{3EA}$$
8. Numerical illustration. For a bar with $EA = 4.2\ \text{MN}$, $L = 2.5\ \text{m}$, $c = 900\ \text{N/m}^{2}$ and $P = 15\ \text{kN}$, the coefficients are $a_1 = 4.241071\times10^{-3}$, $a_2 = 0$, $a_3 = -3.571429\times10^{-5}\ \text{m}^{-2}$, and the tip displacement is $u(L) = 10.0446\ \text{mm}$. The distributed load contributes $cL^{3}/(3EA) = 1.1161\ \text{mm}$ of that, the tip force the remaining $8.9286\ \text{mm}$.