# converting $(a(x)u(x))_{xx}=f(x)$ to weak form in FEM

Suppose I have an ODE of the form

$$u(x) - [a(x)u(x)]_{xx}=f(x)$$

and I want to solve it with FEM, how do I convert it to the weak form? Using integration by part, the middle term becomes

$$(v,-(au)'')=(v',(au)')=(v',a'u+au')=(v',a'u)+(v',au')$$

where $$v=v(x)$$ is the test function. If we let $$v=\phi_i(x)$$, and $$u(x)=\sum_j u_j \phi_j(x)$$, then we end up with something like

$$\sum_j \left[ (\phi_i,\phi_j) + (\phi_i',a'\phi_j) + (\phi_i',a\phi_j')\right] u_j = (\phi_i,f)$$

But I'm not sure if there are any subtleties I have missed. Also now looks like the matrix $$A_{ij}$$ is no longer symmetric, which is different from other more standard problems I have come across like the usual poisson equation $$-(au')'=f$$. So I would like to know how to properly convert equations like this to weak form / solving with FEM.

This problem is often referred to as putting a differential equation in 'divergence form' as $$b(x)u(x)-(c(x)u'(x))' = g(x)$$ since this is a self-adjoint (think symmetric) operator: $$(v,-(cu')') = (v',cu') = (cv',u') = (-(cv')',u)$$ for all sufficiently smooth $$u, v$$. Therefore this form is preferred whenever the weak formulation is used, such as with finite elements or when using an energy method. Your given problem $$u(x) - (a(x)u(x))'' = f(x)$$ is not in divergence form and the resulting operator is not self-adjoint, so we cannot expect the corresponding matrix $$A_{ij}$$ to be symmetric. Sometimes an equation not in divergence form can be converted to one that is, but unfortunately the procedure does not work in this case unless $$a$$ is constant. To see why, we calculate: $$g = bu-(cu')' = bu - c'u' - cu''$$ $$f = u - (au)'' = (1-a'')u - 2a'u' - au''$$ and so equating coefficients of the derivatives of $$u$$ we must have both that $$a = c$$ and $$2a' = c'$$, which is not possible unless $$a' = c' = 0$$ and so $$a$$ is constant. Of course if this were the case then your equation would already be in divergence form by pulling $$a$$ outside of the $$x$$ derivatives.
However, the finite element method still works in the non-symmetric case as long as a weak solution is guaranteed to exist. One can show that this is true for this particular problem using the Lax-Milgram theorem but I'll omit the proof unless specifically asked. But the procedure works exactly the same as usual, just the matrix $$A_{ij}$$ to be inverted is no longer symmetric.