# An Improved Galerkin Framework for Solving Unsteady High-Reynolds Navier-Stokes Equations
## 3 Numerical Methods
### Variational Multiscale Stabilization
A large residual or a large jump term indicates a region where the numerical solution is a poor approximation of the true solution, signaling the need for mesh refinement:  
\[\mathbf{R}_{M}:= \left(\frac{\mathbf{u}_{h}^{n+1}-\mathbf{u}_{h}^{n}}{\Delta t}+ \theta((\mathbf{u}_{h}^{n+1}\cdot\
abla)\mathbf{u}_{h}^{n+1}-\
u\Delta \mathbf{u}_{h}^{n+1}+\
abla p_{h}^{n+1})+(1-\theta)((\mathbf{u}_{h}^{n}\cdot \
abla)\mathbf{u}_{h}^{n}-\right.\] \[\left.\
u\Delta\mathbf{u}_{h}^{n})+\
abla p_{h}^{n}\right)- \mathbf{f}^{n+\theta}+\mathbf{S}_{M}^{n+\theta},\]  
\[\mathbf{R}_{C}:=\
abla\cdot\mathbf{u}_{h}^{n+\theta}+\mathbf{S}_{C}^{n+\theta},\]  
\[a(\mathbf{u},\mathbf{v})=\int_{\Omega}\frac{\partial\mathbf{u}}{\partial t} \cdot\mathbf{v}\,d\Omega+\int_{\Omega}((\mathbf{u}\cdot\
abla)\mathbf{u}) \cdot\mathbf{v}\,d\Omega+\
u\int_{\Omega}\
abla\mathbf{u}\cdot\
abla\mathbf{ v}\,d\Omega,\]  
\[b(\mathbf{u},\mathbf{q})=\int_{\Omega}\mathbf{q}\cdot\
abla\cdot\mathbf{u}\,d\Omega.\]  
It is crucial to note that in the integrated framework, the residuals \(\mathbf{R}_{M}\) and \(\mathbf{R}_{C}\) are defined with respect to the VMS-stabilized equations for the large scales (as implied in (8)), not the original Navier-Stokes equations. This is the cornerstone of the deep coupling between the stabilization and adaptation strategies. The VMS method introduces subgrid-scale models that modify the discrete governing equations. Consequently, a meaningful error estimator must assess the accuracy of the numerical solution \((\mathbf{u}_{h},p_{h})\) against this modified system. By doing so, the AMR strategy is directly informed by the performance of VMS, concentrating mesh refinement in regions where the subgrid-scale model is most active or where the stabilized large-scale solution still exhibits significant residuals. This self-consistent approach ensures that computational effort is focused precisely where it is needed to improve the quality of the stabilized solution. To estimate the error, choose the test functions \(\mathbf{v}=\mathbf{e_{u}}\) and \(\mathbf{q}=\mathbf{e_{p}}\), and substitute them into the error equation, there is  
\[\|\mathbf{e_{u}}\|_{H^{1}}^{2}+\|\mathbf{e_{p}}\|_{L^{2}}^{2}=\int_{\Omega} \mathbf{R}_{M}\cdot\mathbf{e_{u}}\,d\Omega+\int_{\Omega}\mathbf{R}_{C} \mathbf{e_{p}}\,d\Omega+\int_{\partial K}[\
u\partial_{n}\mathbf{u}_{h}-p_{h} \mathbf{n}]\cdot\mathbf{e_{u}}\,dS, \tag{11}\]  
On each cell \(K\), by the Cauchy-Schwarz inequality, one has  
\[\left|\int_{K}\mathbf{R}_{M}\cdot\mathbf{e_{u}}\,dK\right|\leq\|\mathbf{R}_{M }\|_{L^{2}(K)}\|\mathbf{e_{u}}\|_{L^{2}(K)},\]  
\[\left|\int_{K}\mathbf{R}_{C}\mathbf{e_{p}}\,dK\right|\leq\|\mathbf{R}_{C}\|_ {L^{2}(K)}\|\mathbf{e_{p}}\|_{L^{2}(K)}.\]  
\[\left|\int_{\partial K}[\
u\partial_{n}\mathbf{u}_{h}-p_{h}\mathbf{n}]\cdot \mathbf{e_{u}}\,dS\right|\leq\|[\
u\partial_{n}\mathbf{u}_{h}-p_{h}\mathbf{n} ]\|_{L^{2}(\partial K)}\|\mathbf{e_{u}}\|_{L^{2}(\partial K)}.\]Next, based on the properties of the interpolation operator, the error decomposition is given by  
\[\mathbf{R}_{M}:= \left(\frac{\mathbf{u}_{h}^{n+1}-\mathbf{u}_{h}^{n}}{\Delta t}+ \theta((\mathbf{u}_{h}^{n+1}\cdot\
abla)\mathbf{u}_{h}^{n+1}-\
u\Delta \mathbf{u}_{h}^{n+1}+\
abla p_{h}^{n+1})+(1-\theta)((\mathbf{u}_{h}^{n}\cdot \
abla)\mathbf{u}_{h}^{n}-\right.\]0  
where the interpolation error is \(\mathbf{e}_{\mathbf{u}}^{I}=\mathbf{u}-\Pi_{h}\mathbf{u}\), and the discretization error is \(\mathbf{e}_{\mathbf{u}}^{h}=\Pi_{h}\mathbf{u}-\mathbf{u}_{h}\). Then, proceed with the interpolation error estimation. Within element \(K\), a Taylor expansion of \(\mathbf{u}\) is performed up to order \(k+1\) as  
\[\mathbf{R}_{M}:= \left(\frac{\mathbf{u}_{h}^{n+1}-\mathbf{u}_{h}^{n}}{\Delta t}+ \theta((\mathbf{u}_{h}^{n+1}\cdot\
abla)\mathbf{u}_{h}^{n+1}-\
u\Delta \mathbf{u}_{h}^{n+1}+\
abla p_{h}^{n+1})+(1-\theta)((\mathbf{u}_{h}^{n}\cdot \
abla)\mathbf{u}_{h}^{n}-\right.\]1  
where \(\mathbf{x}_{c}\) is the cell center, \(\xi\) is some intermediate point, and \(D\) is the derivative operator.