From 3cb4b7054b3849dfe5369309116f1c69f40239b5 Mon Sep 17 00:00:00 2001
From: David Bold <dave@ipp.mpg.de>
Date: Mon, 15 Jan 2024 13:34:13 +0100
Subject: [PATCH] handle maxfac_inner = 1 differently

If we want them equally spaced we can find a solution much simpler.
---
 zoidberg/poloidal_grid.py | 34 +++++++++++++++++++++++++++++-----
 1 file changed, 29 insertions(+), 5 deletions(-)

diff --git a/zoidberg/poloidal_grid.py b/zoidberg/poloidal_grid.py
index 5fe4633..c87da57 100644
--- a/zoidberg/poloidal_grid.py
+++ b/zoidberg/poloidal_grid.py
@@ -683,7 +683,6 @@ def grid_elliptic(
     https://en.wikipedia.org/wiki/Principles_of_grid_generation
 
     """
-
     if nx_outer:
         assert nx_outer > 0
         nx -= nx_outer
@@ -749,18 +748,43 @@ def dx(x, sign):
         x0 = x.copy()
         for i in range(steps):
             x += fac * (dx(x, 1) + dx(x, -1))
+
+        extra = 0
         if maxfac_inner:
 
             def getfac(x):
                 d = dx(x, -1)
                 return np.max(d) / np.min(d)
 
-            print(f"starting fac {getfac(x)}")
-            while getfac(x) > maxfac_inner:
+            while not np.all(dx(x, -1) > 0):
                 x += fac * (dx(x, 1) + dx(x, -1))
-            print(f"finished fac {getfac(x)}")
+                extra += 1
+
+            if maxfac_inner == 1:
+                x1 = x.copy()
+                for i in range(1, len(x0)):
+                    if x1[i - 1] > x1[i]:
+                        x1[i] += 2 * np.pi
+                x2 = (
+                    np.mean(x1)
+                    + np.linspace(-np.pi, np.pi, len(x0), endpoint=False)
+                    + np.pi / len(x0)
+                )
+                gf = getfac(x0)
+                ratio = (maxfac_inner - 1) / (gf - 1)
+                x = x2 * (1 - ratio) + np.array(x1) * (ratio)
+                x %= 2 * np.pi
+                # print(getfac(x), ratio, getfac(x2), getfac(x1))
+            else:
+                # print(f"starting fac {getfac(x)}")
+                while getfac(x) > maxfac_inner:
+                    extra += 1
+                    # if extra % 100 == 0:
+                    #    print(extra, fac, getfac(x), 0.1 / (getfac(x) - 1))
+                    #    fac = max(0.05 / (getfac(x) - 1), 0.1)
+                    x += fac * (dx(x, 1) + dx(x, -1))
+                # print(f"finished fac {getfac(x)}")
         else:
-            extra = 0
             while not np.all(dx(x, -1) > 0):
                 x += fac * (dx(x, 1) + dx(x, -1))
                 extra += 1