When polynomials of degree \(k\ge 1\)  is used in the scheme,
the global linear system, which is symmetric and indefinite,
consists of degrees of freedom (DOFs) for
the normal component of velocity (of polynomial degree \(k\) ) on the mesh
skeleton (facets), the tangential hybrid velocity (of polynomial degree \(k-1\) )
on the mesh skeleton, and one pressure DOF per element on the
mesh. The linear system problem is then solved via a preconditioned MinRes
method [1]} with a block diagonal preconditioner which is of similar form as the uniform
preconditioner studied in Olshanskii et al. [2]}
for a generalized Stokes interface problem. We further use an auxiliary
space preconditioner of Xu [3]}
with algebraic multigrid (AMG) for the velocity block and a hypre AMG preconditioner
for the pressure block to arrive at the final block AMG preconditioner.
This preconditioner is numerically verified to be robust with respect to
mesh size, time step size, and material parameters.
