1 files changed, 9 insertions, 8 deletions

diff --git a/topology.cpp b/topology.cpp
index 1d46649..129957e 100644
--- a/topology.cpp
+++ b/topology.cpp
@@ -103,8 +103,8 @@ public:
     return {∂L, ∂∂L};
   }
 
-  Vector newtonMethod(const Vector& v0, Real γ = 1, Real ε = 1e-12, unsigned maxSteps = 1e2) const {
-    Vector v = v0;
+  Vector newtonMethod(const Vector& v₀, Real γ = 1, Real ε = 1e-12, unsigned maxSteps = 1e3) const {
+    Vector v = v₀;
     unsigned steps = 0;
 
     Vector ∂L;
@@ -112,6 +112,7 @@ public:
 
     while (std::tie(∂L, ∂∂L) = ∂L_∂∂L(v), ∂L.squaredNorm() > ε) {
       if (v.tail(M + 1).squaredNorm() > 1 / ε || steps > N * maxSteps) {
+        std::cerr << "Quit Newton method algorithm after " << steps << "steps" << std::endl;
         return Vector::Zero(N + M + 1);
       }
 
@@ -124,18 +125,18 @@ public:
     return v;
   }
 
-  Vector randomStationaryPoint(Rng& r, Real γ = 1, Real ε = 1e-12, unsigned maxTries = 1e1) const {
+  Vector randomStationaryPoint(Rng& r, Real γ = 1, Real ε = 1e-12, unsigned maxTries = 1e2) const {
     Vector v = Vector::Zero(N + M + 1);
     unsigned tries = 0;
 
     while (v.squaredNorm() < ε && tries < maxTries) {
-      Vector v0(N + M + 1);
-      for (Real& v0ᵢ : v0) {
-        v0ᵢ = r.variate<Real, std::normal_distribution>();
+      Vector v₀(N + M + 1);
+      for (Real& v₀ᵢ : v₀) {
+        v₀ᵢ = r.variate<Real, std::normal_distribution>();
       }
-      v0.head(N) = normalize(v0.head(N));
+      v₀.head(N) = normalize(v₀.head(N));
 
-      v = newtonMethod(v0, γ, ε);
+      v = newtonMethod(v₀, γ, ε);
       tries++;
     }