Fix MappingQ1

diff --git a/source/fe/mapping_q_generic.cc b/source/fe/mapping_q_generic.cc
index 9f5eec7b3d89dcd1dcc1826ad48e1948cc128e95..4696483c380c04027b5cdf37ac0636f58475c39e 100644 (file)
--- a/source/fe/mapping_q_generic.cc
+++ b/source/fe/mapping_q_generic.cc
@@ -106,8 +106,10 @@ namespace internal
             eta1 = -c/b;
             eta2 = -c/b;
           }
-        // special case #2: if c is very small:
-        else if (std::abs(c/b) < 1e-12)
+        // special case #2: if c is very small or the square root of the
+        // discriminant is nearly b.
+        else if (std::abs(c) < 1e-12*std::abs(b)
+                 || std::abs(std::sqrt(discriminant) - b) <= 1e-14*std::abs(b))
           {
             eta1 = (-b - std::sqrt(discriminant)) / (2*a);
             eta2 = (-b + std::sqrt(discriminant)) / (2*a);