Merge pull request #10 from wellington36/v1.3.0

V1.3.0

.gitignore

@@ -0,0 +1,2 @@
+# python log
+/__pycache__/

README.md

@@ -3,27 +3,27 @@
 This repository contains implementations of the following series transformations:
 
 * [Aitken's transformation (or delta-squared process)](https://en.wikipedia.org/wiki/Aitken%27s_delta-squared_process):
-  - One step: O($n$)
-  - Exhausted: O($n^2$)
+  - Transform: O($n$)
+  - To find n: O($2n\log n$)
 
 * [Richardson's transformation (modify, with given p)](https://en.wikipedia.org/wiki/Richardson_extrapolation):
-  - One step: O($\log n$)
-  - Exhausted: O($(\log n)^2$)
+  - Transform: O($\log n$)
+  - To find n: O($2(\log n)^2$)
 
 * [Epsilon transformation](https://www.sciencedirect.com/science/article/pii/S0377042700003551):
-  - One step: O($n$)
-  - Exhausted: O($n^2$)
+  - Transform: O($n$)
+  - To find n: O($2n\log n$)
 
 * [G transformation](https://epubs.siam.org/doi/abs/10.1137/0704032?journalCode=sjnaam):
-  - One step: O($n$)
-  - Exhausted: O($n^2$)
+  - Transform: O($n$)
+  - To find n: O($2n\log n$)
 
 ## Usage
 
-In `acceleration.py` we have functions that receive a list of the values of the series to be accelerated along with the number of steps (this being a positive integer or -1 if you want to do the process until you can't do it anymore). Returning the new series in the form of a list. For example:
+In `acceleration.py` we have the transformations implemented above, and for use we have the `acceleration` function that receives a series (in the form of a function $s: \mathbb{N} \to \mathbb{R}^n$ returning the first n elements of that series), the chosen transformation ("Aitken_tranform", "Richardson_transform", "Epsilon_transform", "G_transform" e "no_transform", the latter being using the initial series without any transformation), the stopping criterion (in case the difference of the last two values of the series are smaller than a given error) and the maximum number of steps the transformation can take. This function finds the smallest n at which, with the transformation, the series error is less than the given error. For example:
 
 ```python
-from acceleration import Aitken_tranform
+from acceleration import *
 
 # a zeta(2) series
 def square_series(n: int) -> np.ndarray:
@@ -36,6 +36,7 @@ def square_series(n: int) -> np.ndarray:
 
     return series
 
-# create a new series with the 100 first terms using Aitken acceleration
-accelerated_series = Aitken_tranform(square_series(100))
+# Creates a new accelerated series with fewer terms than the original and 
+# such that the difference of the last two terms is less than the error=1e-5:
+accelerated_series = acceleration(square_series, Aitken_tranform, error=1e-5, max_steps=2)
 ```

acceleration.py

@@ -2,9 +2,11 @@
 import numpy as np
 from configuration import *
 
-def Aitken_tranform(items: np.ndarray, steps=-1) -> np.ndarray:
-    if steps == -1:
-        steps = int((len(items) - 1) / 2)
+def no_transform(items: np.ndarray, max_steps=10) -> np.ndarray:
+    return items
+
+def Aitken_tranform(items: np.ndarray, max_steps=10) -> np.ndarray:
+    steps = min(int((len(items) - 1) / 2), max_steps)
 
     for _ in range(steps):
         acel = np.zeros(len(items) - 2, dtype=DT)
@@ -17,13 +19,12 @@ def Aitken_tranform(items: np.ndarray, steps=-1) -> np.ndarray:
 
         if len(acel) < 3:
             return acel
-        
+
     return acel
 
-def Richardson_transform(items: np.ndarray, p=1, steps=-1) -> np.ndarray:
+def Richardson_transform(items: np.ndarray, p=1, max_steps=10) -> np.ndarray:
     """Receive a p that represents the power of the Richardson transform"""
-    if steps == -1:
-        steps = int(math.log2(len(items))) - 1
+    steps = min(int(math.log2(len(items))) - 1, max_steps)
 
     for _ in range(steps):
         acel = np.zeros(int(len(items)/2), dtype=DT)
@@ -37,13 +38,12 @@ def Richardson_transform(items: np.ndarray, p=1, steps=-1) -> np.ndarray:
 
     return acel
 
-def Epsilon_transfom(items: np.ndarray, steps=-1) -> np.ndarray:
+def Epsilon_transform(items: np.ndarray, max_steps=10) -> np.ndarray:
     # Initial values
     aux = np.zeros(len(items)+1, dtype=DT)
     acel = items
-
-    if steps == -1:
-        steps = int(len(items) / 2) - 1
+
+    steps = min(int(len(items) / 2) - 1, max_steps)
 
     for _ in range(steps):
         for i in range(0, len(aux) - 3):
@@ -53,11 +53,10 @@ def Epsilon_transfom(items: np.ndarray, steps=-1) -> np.ndarray:
         for i in range(0, len(acel) - 3):
             acel[i] = acel[i+1] + 1/(aux[i+1] - aux[i])
         acel = acel[:-2]
-
 
     return acel
 
-def G_transform(items: np.ndarray, steps=-1) -> np.ndarray:
+def G_transform(items: np.ndarray, max_steps=10) -> np.ndarray:
     # Initial values
     aux1 = np.ones(len(items) + 1, dtype=DT)
     aux2 = np.zeros(len(items), dtype=DT)
@@ -68,8 +67,7 @@ def G_transform(items: np.ndarray, steps=-1) -> np.ndarray:
 
     acel = items
 
-    if steps == -1:
-        steps = len(items) - 1
+    steps = min(len(items) - 1, max_steps)
 
     for _ in range(steps):
         for i in range(len(aux1) - 2):
@@ -87,6 +85,32 @@ def G_transform(items: np.ndarray, steps=-1) -> np.ndarray:
     return acel
 
 
+def acceleration(series, transform, error=1e-5, max_steps=10) -> np.ndarray:
+    n0 = 10
+    acel = transform(series(n0))
+    i = -1  # trash
+
+    while abs(acel[-1] - acel[-2]) > error: # check error
+        i = i + 1
+        n = n0 + 2**i
+        acel = transform(series(n), max_steps=max_steps)
+
+    n0 = n0 + 2**(i-1)
+
+    while (n > n0):
+        acel = transform(series(int((n+n0)/2)), max_steps=max_steps)
+
+        if abs(acel[-1] - acel[-2]) > error:    # check error
+            n0 = int((n+n0)/2 + 1)
+        else:
+            n = int((n+n0)/2)
+
+    acel = transform(series(n), max_steps=max_steps)
+
+    #print(n)
+    return acel
+
+
 if __name__ == "__main__":
 
     print("Hello World")

configuration.py

@@ -1,4 +1,7 @@
 import numpy as np
 
 # SET DATA TYPE
-DT = np.dtype('float64') # numpy precision
+DT = np.dtype('float64') # numpy precision
+
+# SET zeta(1.2) VALUE APPROXIMATION (~308)
+z = 5.5915824411777507765365631934231432776299032418023310994737282500248979068026510779911033208016101447710090788070826059927576010402832256754907368859609450565260645425520759912470219931989425344148176937499288655408751303431761462302654347948782886616919938182327581584549962638813992208567921853722195176401

test.py

@@ -1,6 +1,6 @@
 import math
 import numpy as np
-from acceleration import Aitken_tranform, Richardson_transform, Epsilon_transfom, G_transform
+from acceleration import *
 from configuration import *
 
 def square_series(n: int) -> np.ndarray:
@@ -9,18 +9,18 @@ def square_series(n: int) -> np.ndarray:
     series[0] = 1.0
 
     for i in range(1, n):
-        series[i] =  series[i-1] + 1/(i+1)**2
+        series[i] =  series[i-1] + 1/(i+1)**(2)
 
     return series
 
-def zeta(n: int, s: float) -> np.ndarray:
-    """Zeta function"""
+def slow_zeta_series(n: int) -> np.ndarray:
+    """Zeta(1.2) series, converges to z (in the other file)"""
     series = np.zeros(n, dtype=DT)
     series[0] = 1.0
 
     for i in range(1, n):
-        series[i] =  series[i-1] + 1/(i+1)**s
-
+        series[i] =  series[i-1] + 1/(i+1)**(1.2)
+    
     return series
 
 def dirichlet_series(n: int) -> np.ndarray:
@@ -35,13 +35,11 @@ def dirichlet_series(n: int) -> np.ndarray:
 
 
 if __name__ == "__main__":
-    STEPS = 2
-    N = 10_000
-    p = 2
-
-    print(math.pi**2 / 6)
-    print(zeta(N, p)[-1])
-    print(Aitken_tranform(zeta(N, p), steps=STEPS)[-1])
-    print(Richardson_transform(zeta(N, p), steps=STEPS)[-1])
-    print(Epsilon_transfom(zeta(N, p), steps=-1)[-1])
-    print(G_transform(zeta(N, p), steps=STEPS)[-1])
+    e = 1e-3
+
+    print(math.pi**2/6)
+    print(acceleration(square_series, no_transform, e)[-1])
+    print(acceleration(square_series, Aitken_tranform, e)[-1])
+    print(acceleration(square_series, Richardson_transform, e)[-1])
+    print(acceleration(square_series, Epsilon_transform, e)[-1])
+    print(acceleration(square_series, G_transform, e)[-1])