Address review comments A bit more comment, but also a more robust protection against a divide by zero, which happens in the colinear double-cusp case, as exercised by tricky-strokes.

diff --git a/shader/flatten.wgsl b/shader/flatten.wgsl
index 0998cc9..6ee46b3 100644
--- a/shader/flatten.wgsl
+++ b/shader/flatten.wgsl

@@ -96,8 +96,8 @@
     if cth0 * cth1 < 0.0 {
         return 2.0;
     }
-    let e0 = (2. / 3.) / (1.0 + cth0);
-    let e1 = (2. / 3.) / (1.0 + cth1);
+    let e0 = (2. / 3.) / max(1.0 + cth0, 1e-9);
+    let e1 = (2. / 3.) / max(1.0 + cth1, 1e-9);
     let s0 = sin(cparams.th0);
     let s1 = sin(cparams.th1);
     let s01 = cth0 * s1 + cth1 * s0;

diff --git a/src/cpu_shader/euler.rs b/src/cpu_shader/euler.rs
index 08c8edd..17c6bc9 100644
--- a/src/cpu_shader/euler.rs
+++ b/src/cpu_shader/euler.rs

@@ -85,10 +85,18 @@
             // Rationale: this happens when fitting a cusp or near-cusp with
             // a near 180 degree u-turn. The actual ES is bounded in that case.
             // Further subdivision won't reduce the angles if actually a cusp.
+            //
+            // A value of 2.0 represents the approximate worst case distance
+            // from an Euler spiral with 0 and pi tangents to the chord. It
+            // is not very critical; doubling the value would result in one more
+            // subdivision in effectively a binary search for the cusp, while too
+            // small a value may result in the actual error exceeding the bound.
             return 2.0;
         }
-        let e0 = (2. / 3.) / (1.0 + cth0);
-        let e1 = (2. / 3.) / (1.0 + cth1);
+        // Protect against divide-by-zero. This happens with a double cusp, so
+        // should in the general case cause subdivisions.
+        let e0 = (2. / 3.) / (1.0 + cth0).max(1e-9);
+        let e1 = (2. / 3.) / (1.0 + cth1).max(1e-9);
         let s0 = self.th0.sin();
         let s1 = self.th1.sin();
         // Note: some other versions take sin of s0 + s1 instead. Those are incorrect.