summaryrefslogtreecommitdiff
path: root/log-fourier_integrator.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'log-fourier_integrator.cpp')
-rw-r--r--log-fourier_integrator.cpp19
1 files changed, 10 insertions, 9 deletions
diff --git a/log-fourier_integrator.cpp b/log-fourier_integrator.cpp
index 7736be8..110b9e4 100644
--- a/log-fourier_integrator.cpp
+++ b/log-fourier_integrator.cpp
@@ -88,9 +88,8 @@ int main(int argc, char* argv[]) {
Real β = 0;
while (β < βₘₐₓ) {
β += Δβ;
- Real μ₁ = μ * 0.99;
+ Real μ₁ = 0;
Real μ₂ = 0;
- μ *= 1.01;
while (true) {
Real ΔC = 100;
Real ΔC₀ = 100;
@@ -147,16 +146,18 @@ int main(int argc, char* argv[]) {
if (std::abs(Cₜ[0] - 1) < ε) {
break;
}
- if (μ₂ < μ₁) {
- /* Cₜ[0] will generally start too big. Until we find a value too small, we
- * scale down proportional too Cₜ[0] */
- if (Cₜ[0] > 1) {
- μ *= sqrt(sqrt(std::abs(Cₜ[0])));
- } else { /* Otherwise, we've found an upper bound for μ */
+ if (μ₁ == 0 || μ₂ == 0) {
+ if (μ₁ == 0 && Cₜ[0] > 1) {
+ /* We found a lower bound */
+ μ₁ = μ;
+ }
+ if (μ₂ == 0 && Cₜ[0] < 1) {
+ /* We found an upper bound */
μ₂ = μ;
}
+ μ *= sqrt(sqrt(std::abs(Cₜ[0])));
} else {
- /* Once μ₂ is set, we can use bisection */
+ /* Once the bounds are set, we can use bisection */
if (Cₜ[0] > 1) {
μ₁ = μ;
μ = μ₁ + (μ₂ - μ₁) / 2;