Local Solvers for High-Order Patch Smoothers via p-Multigrid

I propose a vertex patch smoother where local problems are solved inexactly by a nested, matrix-free p-multigrid, creating a multigrid-within-multigrid framework. A single iteration of the local solver can be evaluated with $\mathcal{O}(p^{d+1})$ operations, and the approach is applicable to non-separable problems on unstructured meshes. Numerical experiments demonstrate limited sensitivity to geometric distortion and high-contrast coefficients. When used in a global geometric multigrid solver, the method achieves robustness with respect to both polynomial degree $p$ and mesh refinement, even on heavily distorted meshes.