  Phase-field fracture models lead to variational problems that can be written
as a coupled variational equality and inequality system. Numerically, such
problems can be treated with Galerkin finite elements and primal-dual active
set methods. Specifically, low-order and high-order finite elements may be
employed, where, for the latter, only few studies exist to date. The most
time-consuming part in the discrete version of the primal-dual active set
(semi-smooth Newton) algorithm consists in the solutions of changing linear
systems arising at each semi-smooth Newton step. We propose a new parallel
matrix-free monolithic multigrid preconditioner for these systems. We provide
two numerical tests, and discuss the performance of the parallel solver
proposed in the paper. Furthermore, we compare our new preconditioner with a
block-AMG preconditioner available in the literature.
