# Bubnov-Galerkin method in 1D: how to handle convective-type nonlinearity?

Consider the BVP: find $u = u(x)$, for $x \in (0,1)$ that satisfies \begin{align} u'' + u u' = f, \\ u'(0) = g_n, u(1) = g_d. \end{align} To derive the weak form for this BVP, we multiply the first equation by a suitably smooth test function $\Phi = \Phi(x)$ and integrate both sides. This leads to \begin{align} - \int_\Omega \Phi ' u' + \int_\Omega \Phi u u' = \int_\Omega \Phi f + \Phi(0) g_n - \Phi(1) u'(1). \end{align} From this equation, we choose our test function space $\mathscr{T}$ to be a subspace of the Sobolev space $H^1(\Omega)$ that eliminate that $\Phi(1) u'(1)$ term, i.e., \begin{align} \mathscr{T} = \{w: w \in H^1(\Omega), w(1) = 0 \}. \end{align} Furthermore, we choose our set of trial functions $\mathscr{S}$ to be a subset of $H^1(\Omega)$ that satisfy the essential boundary condition, i.e., \begin{align} \mathscr{S} = \{v:v \in H^1(\Omega), v(1) = g_d \}. \end{align}

We then consider the WP associated with the original BVP: find $u \in \mathscr{S}$ such that \begin{align} -\int_\Omega \Phi'u'+\int_\Omega \Phi u u' = \int_\Omega \Phi f + \Phi(0)g_n, \qquad \forall \Phi \in \mathscr{T}. \end{align}

For Galerkin methods, we fix a positive integer $n$ and determine an n-dimensional subspace of $\mathscr{T}$ by specifying a basis, i.e., \begin{align} \mathscr{T}_n = \text{span}\{\phi_1, ..., \phi_n\} \subset\mathscr{T}. \end{align}

For Bubnov-Galerkin, we determine an n-dimensional subset of $\mathscr{S}$ by taking $\mathscr{S}_n = \mathscr{T}_n \cup \{g(x)\}$, where $g(x) = g_d$ handles the essential boundary. Our finite element solution' can then be written $u_h(x) = \sum_{j = 1}^n u_j \phi_j(x) + g(x)$.

Substituting $u_h$ into the weak equation should lead to a nonlinear system of equations ($\forall \Phi \in \mathscr{T}$ becomes $\forall \phi_i \in \mathscr{T}_n$): \begin{align} -\int \phi_i \left( \sum_{j = 1}^n u_j \phi_j + g\right)' + \int \phi_i \left( \sum_{j = 1}^n u_j \phi_j + g \right) \left( \sum_{j = 1}^n u_j \phi_j + g \right)' = \int \phi_i f + \phi_i (0). \end{align} The first term on the left-hand-side can be split up into \begin{align} - \sum_{j = 1}^n \int u_j \phi_i \phi_j' - \int \phi_i g', \end{align} which is just a bilinear form (matrix) and linear form (vector). The entire right hand side consists of just linear forms (vectors).

What do I do with the term in the middle?

If I distribute or foil'' the middle term, I get \begin{align} \int \phi_i \left(\sum_{j = 1}^n u_j \phi_j \right)\left(\sum_{j = 1}^n u_j \phi_j' \right) + \int \phi_i \left(\sum_{j = 1}^n u_j \phi_j g' \right) + \int \phi_i \left(g \sum_{j = 1}^n u_j \phi_j' \right) + \int \phi_i g g'. \end{align} Here the second and third terms look like bilinear forms (matrices), the fourth term looks like a linear form (vector), but what about the first term, the "trilinear form"? It seems to me that this first term would end up like $u^T A u$ since the unknown coefficients $u_j$ appear twice, but I am not too sure.

Also, my motivation for this problem is to understand similar looking nonlinearities are handled practically, because they seem to appear a lot (e.g., in Navier-Stokes).

## 1 Answer

You're on the right track -- the term $$\int \phi_i u u'$$ can be expanded to $$\int \phi_i \left(\sum_j u_j \phi_j\right) \left(\sum_k u_k \phi_k'\right).$$ (Note how I use a different "silent" index $j$ and $k$ in the two expansions.) This leads to a trilinear form $$B_{ijk} = \int \phi_i \phi_j \phi_k'$$ and the term in the variational formulation would be $(BU)U$ with summation over the last index of the tensor and the only index of the vector of coefficients $U$.

Of course, all of this does not help you in practice. We do not know how to solve general quadratic problems in multiple variables other than by resorting to linearization and solving a sequence of linear problems. In other words, we cannot directly solve the problem you have, but only a sequence of problems that are linear in a solution variable. These linear problems must all have bilinear forms of the traditional form.

• Do you know resources for solving such problems in practice? I have seen that one can use Newton's method on $g(u) = u'' + u u' - f$, with initialization $u_0$ satisfies boundary-condition. Then you get $u_{i + 1} = u_i - [\nabla g(u_i)]^{-1} [g(u_i)] = u_i - S$. To determine $S$, we try to solve $[\nabla g(u_i)](s) = -[g(u_i)]$. This is equivalent to solving a linear problem $S'' + S' u_i + u_i' S = -u_i'' - u_i u_i' + f$ in my case. Then this can be solved "standardly". Does this work well enough, or are there better ways? – amarney Sep 29 '17 at 23:08
• Yes, that's the way to do it if you have a good enough initial guess. You can find a number of ways to solve nonlinear problems discussed in lectures 31.5 and following here: math.colostate.edu/~bangerth/videos.html – Wolfgang Bangerth Sep 30 '17 at 14:33
• @WolfgangBangerth I asked a related question here, following this idea and using linear finite elements in 1D scicomp.stackexchange.com/questions/35967/… – Vefhug Sep 23 '20 at 13:32