This list of solvers that have worked fairly consistently, try them out locally and if the recovery solver is not activated, you should be safe for longer transient runs:

% biconjugate gradient with block Jacobi preconditioner [Thiago's favorite]
md.toolkits.DefaultAnalysis=bcgslbjacobioptions();

% gnres with block jacobi preconditioner [Martin's favorite]
md.toolkits.DefaultAnalysis=gmresbjacobioptions();
md.toolkits.DefaultAnalysis.ksp_max_it=1;

% ASM, used for mass conservation but may work for other things
md.toolkits.DefaultAnalysis=asmoptions();

% An iterative solver may work for one analysis and not another, you can assign a solver to each, like this:
md.toolkits.StressbalanceAnalysis = gmresbjacobioptions();
md.toolkits.DefaultAnalysis=issmmumpssolver();

Make sure to always keep MUMPS as the recovery solver (it is slow but extremely robust):

%Keep recovery mode as MUMPS
md.toolkits.RecoveryAnalysis=issmmumpssolver();