I have updated the volcano / dike example on the PyLith wiki. In smoothing the mesh using CUBIT 15.3 (and presumably current versions of Trelis as well), the mesh topology became incorrect resulting in NaN for the norm of the residual when running PyLith.
Simply multiple smoothing passes with a single pass creates a good quality mesh with correct topology. I also updated the field split solver settings in
pylithapp.cfg for PyLith v2.2.1.
Thanks to Josh Jones for reporting this issue.