0.0.2.4

Optimized GLpoint implementation

Author: iparsw

CHANGELOG.md

@@ -6,10 +6,31 @@
 
 bench_1 Result:
 
-| time |  RL         | GL          |
-| ---- | ----------- | ----------- |
-| 1+e2 | 0.0482 ms   | *           |
-| 1+e3 | 5.6009 ms   | 0.9822 ms   |
-| 1+e4 | 605.5952 ms | 5.2982 ms   |
-| 1+e5 | *           | 104.1740 ms |
-| 1+e6 | *           | 667.7444 ms |
+| time | GL          |  RL         | GLpoint     | RLpoint    |
+| ---- | ----------- | ----------- | ----------- | ---------- |
+| 1+e2 | *           | 0.0482 ms   | *           | 0.0099 ms  |
+| 1+e3 | 0.9822 ms   | 5.6009 ms   | 0.0083 ms   | 0.0899 ms  |
+| 1+e4 | 5.2982 ms   | 605.5952 ms | 0.0721 ms   | 1.0599 ms  |
+| 1+e5 | 104.1740 ms | *           | 1.5568 ms   | 9.6303 ms  |
+| 1+e6 | 667.7444 ms | *           | 10.0608 ms  | 93.9401 ms |
+
+
+## 0.0.2.4: Optimized GLpoint implementation
+
+- GLpoint
+  - On-the-fly Coefficient Calculation
+  - Efficient Memory Access
+  - Precomputed Constants
+- GLcoeffs optimization
+
+
+bench_1 Result and comparison to original differint package (1.0.0):
+
+| time |  0.0.2.3    |  0.0.2.4    | differint   |
+| ---- | ----------- | ----------- | ----------- |
+| 1+e2 | *           | *           | *           |
+| 1+e3 | 0.0083 ms   | 0.0099 ms   | 0.4213 ms   |
+| 1+e4 | 0.0721 ms   | 0.0243 ms   | 3.8376 ms   |
+| 1+e5 | 1.5568 ms   | 0.1885 ms   | 36.3030 ms  |
+| 1+e6 | 10.0608 ms  | 3.9995 ms   | 385.9205 ms |
+| 1+e7 | *           | 41.2 ms     | *           |

README.MD

@@ -1,8 +1,10 @@
+- **Only Windows Support for now!**
+
 # differintC
 
 `differintC` is a high-performance C++ library with Python bindings for computing fractional differintegrals (derivatives and integrals of arbitrary real order) using numerical methods.
 
-This package implements optimized versions of the Riemann–Liouville and Grünwald–Letnikov (GL) fractional differintegral operators, inspired by the original [DifferInt project](https://placeholder-link-to-original-differint-project).
+This package implements optimized versions of the Riemann–Liouville and Grünwald–Letnikov (GL) fractional differintegral operators, inspired by the original [DifferInt project](https://github.com/differint/differint).
 
 > ⚙️ Built with modern C++17 and exposed to Python via `pybind11`, this library is significantly faster than pure-Python equivalents, especially for large arrays and high-precision needs.
 
@@ -14,13 +16,6 @@ This package implements optimized versions of the Riemann–Liouville and Grünw
 pip install differintC
 ````
 
-To build from source:
-
-```bash
-git clone https://github.com/your-username/differintC-project.git
-cd differintC-project
-pip install .
-```
 
 ---
 

differintC/_version.py

@@ -1 +1 @@
-__version__ = "0.0.2.3"  # <-- ANY PEP 440 version string allowed
+__version__ = "0.0.2.4"  # <-- ANY PEP 440 version string allowed

src/differint.cpp

@@ -131,16 +131,17 @@ std::vector<T> RL(T alpha,
 // Compute GL coefficients up to order n
 template <typename T>
 std::vector<T> GLcoeffs(T alpha, std::size_t n) {
-    // b[0] = 1
-    std::vector<T> b(n + 1, T(1));
+    std::vector<T> b;
+    b.reserve(n + 1);  // Preallocate memory
+    b.push_back(T(1));  // b0 = 1
+
     for (std::size_t j = 1; j <= n; ++j) {
-        b[j] = b[j - 1] * (static_cast<T>(-alpha) + static_cast<T>(j - 1))
-               / static_cast<T>(j);
+        T j_t = static_cast<T>(j);
+        b.push_back(b.back() * (j_t - 1 - alpha) / j_t);
     }
     return b;
 }
 
-
 // GL at a single point (endpoint)
 template <typename T>
 T GLpoint(T alpha,
@@ -158,15 +159,31 @@ T GLpoint(T alpha,
     if (f_vals.size() != num_points) {
         throw std::invalid_argument("f_vals size must equal num_points");
     }
-    T step = (domain_end - domain_start) / static_cast<T>(num_points - 1);
-    // Compute coefficients for k = num_points - 1
-    std::size_t k = num_points - 1;
-    auto b = GLcoeffs(alpha, k);
-    T acc = T(0);
-    for (std::size_t j = 0; j <= k; ++j) {
-        acc += b[j] * f_vals[k - j];
+
+    const T step = (domain_end - domain_start) / static_cast<T>(num_points - 1);
+    const T step_power = std::pow(step, -alpha);
+    const std::size_t k = num_points - 1;
+
+    T acc = 0;
+    T c_val = 1.0;  // Initialize to b0 = 1
+    const T* f_ptr = f_vals.data() + k;  // Pointer to last element (f_vals[k])
+
+    for (std::size_t j_index = 0; j_index <= k; ++j_index) {
+        // Add current term: c_val * f_vals[k - j_index]
+        acc += c_val * (*f_ptr);
+
+        // Prepare for next iteration (skip for last element)
+        if (j_index < k) {
+            // Update coefficient using recurrence relation:
+            // c_{j+1} = c_j * ( -alpha + j_index ) / (j_index + 1)
+            c_val *= ( -alpha + static_cast<T>(j_index) ) / static_cast<T>(j_index + 1);
+        }
+
+        // Move pointer backward through f_vals
+        --f_ptr;
     }
-    return std::pow(step, -alpha) * acc;
+
+    return step_power * acc;
 }
 
 // GL over entire grid