update to prints, more numba, todos and further bug discoveries

Author: aetherspritee

CCA.py

@@ -20,6 +20,9 @@ def CCA_subcluster(R: np.ndarray, N: int, DF: float, kf: float,iter: int,N_subcl
 
     PCA_OK, data, n_clusters, i_orden = PCA_subcluster(N, N_subcluster, R, DF, kf, tolerance)
 
+    print("==============================")
+    print("         PCA DONE!            ")
+    print("==============================")
     if not PCA_OK:
         return CCA_OK, PCA_OK
 

PCA.py

@@ -1,4 +1,5 @@
 import numpy as np
+import numba
 import time
 
 def PCA(number: int, mass: np.ndarray, r: np.ndarray, Df: float, kf: float, tolerance: float) -> tuple[bool, np.ndarray]:
@@ -89,6 +90,7 @@ def PCA(number: int, mass: np.ndarray, r: np.ndarray, Df: float, kf: float, tole
             m1 = m3
             rg1 = (np.exp(np.sum(np.log(r))/np.log(r).size))*(np.power(n1/kf, 1/Df))
             k = k + 1
+            # FIXME: this gets stuck for some reason (no apparent connection to the arcos issue)
 
     data_new = np.zeros((number,4))
     for i in range(number):
@@ -179,6 +181,7 @@ def gamma_calc(rg1: float,rg2: float,rg3: float,n1: int,n2: int,n3: int) -> tupl
 
     return gamma_ok, gamma_pc
 
+@numba.njit()
 def random_list_selection(gamma_ok: bool, gamma_pc: float,X: np.ndarray, Y: np.ndarray, Z: np.ndarray,R: np.ndarray, n1: int, x_cm: float, y_cm: float, z_cm: float) -> tuple[np.ndarray, float]:
     candidates = np.zeros((n1))
     rmax = 0.0
@@ -226,6 +229,7 @@ def search_monomer_candidates(R: np.ndarray, M: np.ndarray, monomer_candidates:
 
     candidates, rmax = random_list_selection(gamma_ok, gamma_pc, X,Y,Z,R,n1,x_cm,y_cm,z_cm)
 
+    # FIXME: this fails once kf > 1.5 LMFAO
     candidates[RS_1-k-1] = 1
     return candidates, rmax
 
@@ -270,6 +274,10 @@ def sticking_process(x: float,y: float,z: float,r: float,r_k: float, x_cm: float
     distance = np.sqrt(np.power(x2-x1,2) + np.power(y2-y1,2) + np.power(z2-z1,2))
 
     alpha = np.arccos((np.power(r1,2) + np.power(distance,2) - np.power(r2,2))/(2*r1*distance))
+    if abs((np.power(r1,2) + np.power(distance,2) - np.power(r2,2))/(2*r1*distance)) > 1:
+        # FIXME: this should never happen!
+        print(f"{(np.power(r1,2) + np.power(distance,2) - np.power(r2,2))/(2*r1*distance) = }")
+        exit(-1)
     r0 = r1*np.sin(alpha)
 
     # AmBdC = (A+B)/C

README.md

@@ -1 +1,12 @@
 # TunableCCA
+Implementation of FracVAL [[[https://www.sciencedirect.com/science/article/pii/S0010465519300323?via%3Dihub][Morán, J. et al. 2019]]] in python.
+
+## TODO
+ - [ ] Fix PCA bugs
+   + [ ] Instability with arcos arguments, tied to monomer size somehow
+   + [ ] low fractal dimensions fail
+   + [ ] high fractal prefactors fail
+   + [ ] yet to be understood infinite loop occurs sometimes
+ - [ ] Allow distribution functions for monomer radii
+ - [ ] Parallelize
+ - [ ] High fractal dimensions are super slow, investigate why

main.py

@@ -6,9 +6,9 @@
 from CCA import CCA_subcluster
 
 # config
-DF = 2.0
+DF = 2.7
 Kf = 1.0
-N = 1024
+N = 64
 R0 = 0.01
 SIGMA = 0
 EXT_CASE = 0