# SUPG Stabilization For streamline upwind/Petrov-Galerkin (SUPG) (see for example {cite:t}john:knobloch:2006,clevenger:heister:2019), we add to the weak form $a(\cdot,\cdot)$ the cell-wise defined weak form {math} a_{\text{SUPG}} (T, \varphi) = \sum_{K \in \mathcal{T}_h} \delta_K \left( \rho C_p \frac{\partial T}{\partial t} - k \triangle T + \mathbf{\beta} \cdot \nabla T - F, \mathbf{\beta} \cdot \nabla \varphi \right)_K,  where $K \in \mathcal{T}_h$ are the cells in the computation, $\delta_K \geq 0$ is a stabilization coefficient defined on each cell, $\mathbf{\beta} = \rho C_p \mathbf{u}$ is the effective advection velocity. The standard literature about SUPG does not contain $\rho C_p$, so it makes sense to include this in the velocity. The first argument in the inner product is the strong form of the residual of PDE, which is tested with the expression $\mathbf{\beta} \cdot \nabla \varphi$ representing the solution in streamline direction. We have to assume $k$ to be constant per cell, as we can not compute the spatial derivatives easily. For the implementation, $\frac{\partial T}{\partial t}$ is replaced by the BDF2 approximation, and its terms from older timesteps and $-F$, are moved to the right-hand side of the PDE. We use the parameter design presented in {cite:t}john:knobloch:2006 for $\delta_K$: {math} \delta_K = \frac{h}{2d\|\mathbf{\beta}\|_{\infty,K}} \left( \coth(Pe)-\frac{1}{Pe} \right)  where the Peclet number is given by {math} Pe = \frac{ h \| \mathbf{\beta} \|_{\infty,K}}{2 d k_\text{max}},  $d$ is the polynomial degree of the temperature or composition element (typically 2), $\coth(x) = (1+\exp(-2x)) / (1-\exp(-2x)),$ and $k_\text{max}=\| k \|_{\infty, K}$ is the maximum conductivity in the cell $K$. If $Pe<1$, the equation is diffusion-dominated and no stabilization is needed, so we set $\delta_K=0$. Care needs to be taken in the definition if $\| \beta \|$ or $k$ become zero: 1. If $k$ is zero, then $Pe=\infty$ and the right part of the product in the definition of $\delta_K$ is equal to one. 2. If $\| \beta \|$ is zero, $Pe < 1$, so we set $\delta_K=0$. 3. If both are zero, no stabilization is needed (the field remains constant).