From a949652cb1aa74ce4aa892392e49de1c32ce32da Mon Sep 17 00:00:00 2001
From: Amaury Bilocq <amaury.bilocq@uliege.be>
Date: Wed, 5 Mar 2025 21:13:12 +0100
Subject: [PATCH] Split TKE contributions
---
pyTurbulence/solver.py | 18 +++++++++++-------
pyTurbulence/spectrum.py | 4 ++--
2 files changed, 13 insertions(+), 9 deletions(-)
diff --git a/pyTurbulence/solver.py b/pyTurbulence/solver.py
index d8cdbb4..8fd62e2 100644
--- a/pyTurbulence/solver.py
+++ b/pyTurbulence/solver.py
@@ -55,15 +55,15 @@ def solve(user_params)->dict:
# Compute additional parameters
mean_density = mean_pressure / mean_temperature
celerity = np.sqrt(gamma * mean_pressure / mean_density)
- Ek = 0.5 * (Mt * celerity) ** 2
- urms = np.sqrt(2. * Ek / 3.)
+ TKE_expected = 0.5 * (Mt * celerity) ** 2
+ urms = np.sqrt(2. * TKE_expected / 3.)
# Generate the solenoidal velocity field
start = time.time()
us, vs, ws, wmax = compute_solenoidal_velocities(spectrum, nbModes, urms,
k0, domain_size,
domain_resolution)
- TKE = np.mean(0.5 * (us.reshape(-1) ** 2 + vs.reshape(-1) ** 2 + ws.reshape(-1) ** 2))
+ TKE_solenoidal = np.mean(0.5 * (us.reshape(-1) ** 2 + vs.reshape(-1) ** 2 + ws.reshape(-1) ** 2))
end = time.time()
solenoidal_time = end - start
@@ -92,6 +92,9 @@ def solve(user_params)->dict:
ud,vd,wd,
incompressible_pressure_fluctuations,
gamma, Mt, params["case"])
+ TKE_total = np.mean(0.5 * (u.reshape(-1) ** 2 + v.reshape(-1) ** 2 + w.reshape(-1) ** 2))
+ TKE_dilatational = TKE_total - TKE_solenoidal
+ TKE = [TKE_solenoidal, TKE_dilatational, TKE_total, TKE_expected]
end = time.time()
thermodynamic_time = end - start
@@ -133,7 +136,6 @@ def solve(user_params)->dict:
results = {
"TKE": TKE,
- "Ek": Ek,
"solenoidal_time": solenoidal_time,
"pressure_time": pressure_time,
"dilatational_time": dilatational_time,
@@ -225,7 +227,7 @@ def _validate_params(params) -> dict:
"Pressure": (0, None), # Pressure must be positive
"Temperature": (0, None), # Temperature must be positive
"k0": (1, None), # Characteristic wavenumber must be >= 1
- "gamma": (0, None), # Heat capacity ratio is typically between 1 and 2
+ "gamma": (0, None) # Heat capacity ratio is typically between 1 and 2
}
validated_params = default_params.copy()
@@ -321,8 +323,10 @@ def _write_log_file(params, results):
f.write(f"{'Plotting the spectrum':<50} {results['spectrum_time']:>20.4f}\n")
f.write(f"{'Saving the data':<50} {results['save_time']:>20.4f}\n")
f.write(f"{'-'*70}\n")
- f.write(f"{'Total TKE':<50} {results['TKE']:>20.4f}\n")
- f.write(f"{'Expected TKE':<50} {results['Ek']:>20.4f}\n")
+ f.write(f"{'Solenoidal TKE':<50} {results['TKE'][0]:>20.4f}\n")
+ f.write(f"{'Dilatational TKE':<50} {results['TKE'][1]:>20.4f}\n")
+ f.write(f"{'Total TKE':<50} {results['TKE'][2]:>20.4f}\n")
+ f.write(f"{'Expected TKE':<50} {results['TKE'][3]:>20.4f}\n")
f.write(f"{'='*70}\n")
# Print log to console
diff --git a/pyTurbulence/spectrum.py b/pyTurbulence/spectrum.py
index a1f2406..313537c 100644
--- a/pyTurbulence/spectrum.py
+++ b/pyTurbulence/spectrum.py
@@ -60,8 +60,8 @@ def compute_tke_spectrum(u : np.ndarray,v : np.ndarray,w : np.ndarray,
# --- Size
n0, n1, n2 = u.shape
- L = [lx, ly, lz]
- oneOverN = 1./(n0*n1*n2)
+ L = [lx, ly, lz]
+ oneOverN = 1./(n0*n1*n2)
# --- FFT
fourier_u = np.fft.fftn(u)*oneOverN
--
GitLab