diff options
author | Jaron Kent-Dobias <jaron@kent-dobias.com> | 2024-08-25 17:08:52 +0200 |
---|---|---|
committer | Jaron Kent-Dobias <jaron@kent-dobias.com> | 2024-08-25 17:08:52 +0200 |
commit | e1ac241062fbe4632e58ec26052d152f6e5aba4b (patch) | |
tree | 35c654cacca8645164eec0852e5ea20cd7b6d164 /topology.cpp | |
parent | 5dca01494df22217e51457377c6d055a6125d647 (diff) | |
download | code-e1ac241062fbe4632e58ec26052d152f6e5aba4b.tar.gz code-e1ac241062fbe4632e58ec26052d152f6e5aba4b.tar.bz2 code-e1ac241062fbe4632e58ec26052d152f6e5aba4b.zip |
Diffstat (limited to 'topology.cpp')
-rw-r--r-- | topology.cpp | 17 |
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++; } |