Remove seed generator and add description of input parameters

example.py

@@ -15,7 +15,6 @@ def main():
         "k0": 12.0,
         "domain_size": (2. * np.pi, 2. * np.pi, 2. * np.pi),
         "domain_resolution": (64,64,64),
-        "seed": 42,
         "gamma": 1.4,
         "case": 1,
         "Spectrum": "PassotPouquet",

pyTurbulence/solver.py

@@ -12,8 +12,20 @@ def solve(user_params)->dict:
 
     Parameters
     ----------
-    user_params : dict
-        Dictionary containing the simulation parameters.
+    params : dict
+        Dictionary containing the simulation parameters. Expected keys and descriptions:
+        
+        - "nbModes" (int): Number of modes per wave (default: 128).
+        - "Mt" (float): Turbulent Mach number (default: 0.2).
+        - "Pressure" (float): Mean pressure in the simulation (default: 1.0).
+        - "Temperature" (float): Mean temperature in the simulation (default: 1.0).
+        - "k0" (int): Peak wavenumber (default: 4).
+        - "domain_size" (list, tuple, np.ndarray): Size of the simulation domain, should be a 3D vector (default: [2π, 2π, 2π]).
+        - "domain_resolution" (list, tuple, np.ndarray): Resolution of the domain, should be a 3D vector (default: [64, 64, 64]).
+        - "Spectrum" (str): Type of spectrum for turbulence generation. Options: "PassotPouquet", "Constant" (default: "PassotPouquet").
+        - "gamma" (float): Heat capacity ratio, typically between 1 and 2 (default: 1.4).
+        - "case" (int): Dimensionalisation of the density and temperature fluctuations (from Ristorcelli & Blaisdell). Options: 1, 2 (default: 1).
+        - "output_dir" (str): Directory to store output results (default: "results").
 
     Returns
     -------
@@ -48,7 +60,7 @@ def solve(user_params)->dict:
 
     # Generate the solenoidal velocity field
     start = time.time()
-    u, v, w, wmax = compute_solenoidal_velocities(spectrum, nbModes, urms, k0, domain_size, domain_resolution, seed=params["seed"])
+    u, v, w, wmax = compute_solenoidal_velocities(spectrum, nbModes, urms, k0, domain_size, domain_resolution)
     TKE = np.mean(0.5 * (u.reshape(-1) ** 2 + v.reshape(-1) ** 2 + w.reshape(-1) ** 2))
     end = time.time()
     solenoidal_time = end - start
@@ -149,7 +161,6 @@ def _validate_params(params) -> dict:
         "domain_resolution": [64, 64, 64],  
         "Spectrum": "PassotPouquet",  
         "gamma": 1.4,  
-        "seed": 42,  
         "case": 1,
         "output_dir": "results"
     }
@@ -165,7 +176,6 @@ def _validate_params(params) -> dict:
         "domain_resolution": (list, tuple, np.ndarray),
         "Spectrum": str,
         "gamma": (int, float),
-        "seed": int,
         "case": int,
         "output_dir": str
     }

pyTurbulence/syntheticTurbulence.py

@@ -24,8 +24,8 @@ def _compute_velocity_field(u_, v_, w_, nx, ny, nz, kx, ky, kz, xc, yc, zc, um,
 # Modified by Amaury Bilocq in 2024
 def compute_solenoidal_velocities(spectrum: str, nbModes: int, urms: float,
                                     k0: float, domain_size: Tuple[float, float, float], 
-                                    domain_resolution: Tuple[int, int, int],
-                                    seed=None) -> Tuple[np.ndarray, np.ndarray, np.ndarray, float]:
+                                    domain_resolution: Tuple[int, int, int]
+                                    ) -> Tuple[np.ndarray, np.ndarray, np.ndarray, float]:
 
     """
     Generate a synthetic turbulence velocity field using the Tony Saad methodology.
@@ -44,8 +44,6 @@ def compute_solenoidal_velocities(spectrum: str, nbModes: int, urms: float,
         Physical size of the domain in each direction (Lx, Ly, Lz).
     domain_resolution : np.ndarray
         Number of grid points in each direction (nx, ny, nz).
-    seed : int, optional
-        Seed for the random number generator.
 
     Returns
     -------
@@ -58,8 +56,6 @@ def compute_solenoidal_velocities(spectrum: str, nbModes: int, urms: float,
     wmax : float
         Maximum wavenumber.
     """
-    if seed is not None:
-        np.random.seed(seed)
 
     Lx, Ly, Lz = domain_size
     nx, ny, nz = domain_resolution

test/test_SolenoidalVelocity.py

@@ -20,17 +20,17 @@ def test_solenoidal_velocity():
     domain_size = (2. * np.pi, 2. * np.pi, 2. * np.pi)
     domain_resolution = (64,64,64)
 
-    u,v,w,wmax = compute_solenoidal_velocities(name,nbModes, urms, k0, domain_size, domain_resolution, seed=32)
+    u,v,w,_ = compute_solenoidal_velocities(name,nbModes, urms, k0, domain_size, domain_resolution)
 
     # Compute TKE from generated field
     TKE_field = np.mean(0.5*(u.reshape(-1)**2 + v.reshape(-1)**2 + w.reshape(-1)**2))
 
     # Generate spectrum from solenoidal field
-    knyquist, wave_numbers, E_k = compute_tke_spectrum(u, v, w, domain_size[0], domain_size[1], domain_size[2])
+    _, wave_numbers, E_k = compute_tke_spectrum(u, v, w, domain_size[0], domain_size[1], domain_size[2])
     idxMax = np.argmax(E_k)
     k0_field = wave_numbers[idxMax]
 
-    assert np.isclose(TKE_field, TKE_spectrum, rtol=0.005), "Test failed: Generated field does not match expected TKE"
+    assert np.isclose(TKE_field, TKE_spectrum, rtol=0.05), "Test failed: Generated field does not match expected TKE"
     assert np.isclose(k0_field, k0_spectrum, rtol=0.01), "Test failed: Generated field does not conserve TKE"