diff --git a/dune/gfe/trustregionmmgbasesolver.hh b/dune/gfe/trustregionmmgbasesolver.hh
index 86776e0fed5ead24733c839b2fcdf04e4e0e7b0a..50c0522ad364d264106ce5d64deaadf095187fd6 100644
--- a/dune/gfe/trustregionmmgbasesolver.hh
+++ b/dune/gfe/trustregionmmgbasesolver.hh
@@ -138,7 +138,7 @@ void TrustRegionMMGBaseSolver<MatrixType, VectorType>::solve()
   matrix_->umv(*x_, tmp);
   double modelDecrease = (modifiedRhs*(*x_)) - 0.5 * ((*x_)*tmp);
 
-  if (modelDecrease < 0)
+  if (std::isnan(modelDecrease) or modelDecrease < 0)
   {
     std::cout << "Model increase: " << -modelDecrease << ", falling back to slower solver" << std::endl;