Extension of Holographic Cosmology to Higher Dimensions and the Dimensional Scale Invariance of Entropic Forces

Extension of Holographic Cosmology to Higher Dimensions and the Dimensional Scale Invariance of Entropic Forces Daisuke SATO1,2* 1*Comprehensive Research Organization for Science and Society, Tsukuba Industry-Academic Collaboration Building, 1601 Kamitakatsu, Tsuchiura City, Ibaraki Prefecture, JAPAN. 2College of Science, Engineering and Technology, University of South Africa, NB Pityina Building Florida, Johannesburg, Gauteng, Republic of South Africa. Corresponding author(s). E-mail(s): [email protected]; ORCID: 0009-0008-3878-4169; Abstract In this study, I present an extension to arbitrary Ddimensions, and through precise dimensional analysis of holographic screens in higher-dimensional Ddimensional spacetime, I clarify that the area scaling A(L, D) = A0LD−2, the information density σscreen (L, D)=σ0/LD−2, the dimensional invariance of the entropic force F=Ts dS dx in arbitrary dimensions, and the scale invariance under rescaling L→λL are mathematically extensible to arbitrary dimensions. Keywords: Cosmology, Gravitational Thermodynamics, Thermodynamics, Gravity, Entropy Growth, Non-equilibrium Structures, Holographic thermodynamics system, 1 Introduction The extension of entropic force theory to higher dimensions provides a pathway to unify string theory, brane theory, and supergravity theory, enabling fundamental theoretical developments in holographic cosmology [4,10,11]. This section presents a rigorous theoretical construction based on dimensional analysis, discusses potential connections to AdS/CFT correspondence, and verifies strict consistency with the 1 holographic principle that entropy S is proportional to area A (and invariant under rescaling) [1,2]. 2 Higher-Dimensional Holographic Screens: Dimensional Analysis 2.1 Basic Geometric Scaling In arbitrary D -dimensional spacetime, holographic screens are defined as (D−1) -dimensional hypersurfaces [1,2], with spatial cross-sections possessing (D−2) dimensions. The scaling relations derived from this geometric structure are Area Scaling Law A(L, D) = A0·LD−2(1) Information Density Scaling σscreen(L, D) = σ0 LD−2(2) where A0 and σ0 are dimension-independent constants [3,23]. 2.2 Dimensional Invariance of Entropic Force The fundamental equation for entropic force is F=Ts dS dx (3) Verification through Dimensional Analysis Temperature dimension [Ts] = K 2 (invariant across arbitrary dimensions) Entropy gradient dimension Dimensionless entropy [S] = 1 (dimensionless through product of information density and area) Gradient dS dx =m−1 Force dimension [F]=[Ts]·dS dx =K·m−1(4) Including the Boltzmann constant [F] = kB·K·m−1=J K·K·m−1=J·m−1=kg ·m·s−2(5) Important Consequence For arbitrary dimension D , through appropriate information density scaling σ∝L−(D−2) , entropic force maintains the physically meaningful force dimension [kg ·m·s−2] [3,8]. 2.3 Consistency with Holographic Principle and Scale Invariance The theoretical framework strictly aligns with the requirement that entropy S is proportional to area A (and invariant under rescaling), as demonstrated through mathematical derivation and dimensional analysis. The proof leverages dimensional analysis [1,2] and specific examples, grounded in the holographic principle’s scale invariance [18]. 3 2.3.1 Theoretical Foundations In holographic theory, entropy satisfies S∝A (proportional to area, independent of volume) [18]. Given the scalings A∝LD−2 and σ∝L−(D−2) , the total entropy is S=σA , which must be constant (invariant) under length rescaling L→λL [1,2]. 2.3.2 Scale Transformation Derivation Under the length scale transformation L→λL •Area transformation A(λL) = λD−2A(L) •Information density transformation σ(λL) = λ−(D−2)σ(L) •Total entropy S(λL) = σ(λL)·A(λL) = λ−(D−2) ·λD−2·S(L) = S(L) Thus, entropy S becomes truly scale-invariant, guaranteeing the physical consistency of entropic force across arbitrary dimensions [3]. More explicitly •Entropy Expression S=L−(D−2) ·LD−2= 1 (6) (constant, independent of L 4 ) [1]. Length Rescaling L→λL S(λL) = (λL)−(D−2) ·(λL)D−2=λ−(D−2) ·λD−2·S(L) = S(L)(7) (invariant) [1]. The ratio S(λL)/S(L) = 1 [2]. Gradient ∂LS= 0 (8) (consequence of invariance) [23]. •Force Dimension F=Ts∂xS , with [Ts] = K and [∂xS] = m−1 , yields [F] = kg ·m·s−2 (consistent) [4]. 2.4 Theoretical Derivation of the Stefan-Boltzmann Law in D-Dimensional Spacetime To provide a rigorous theoretical foundation for the thermodynamic scaling relation u∝T12 introduced in the main text for D= 12 dimensions, we present here the general derivation of the Stefan-Boltzmann law in arbitrary D-dimensional spacetime (comprising 1 time dimension and D−1spatial dimensions). 2.4.1 Generalization of the Planck Distribution For photons with energy E=ℏω, the Bose-Einstein distribution is given by: n(ω) = 1 eℏω/(kBT)−1(9) 2.4.2 Density of States in (D−1)-Dimensional Space In (D−1)-dimensional spatial manifolds, the density of states in momentum space is determined by the volume of a hyperspherical shell: g(k)∝kD−2dk (10) Using the dispersion relation ω=ck for massless photons, we obtain: g(ω)∝ωD−2dω (11) 5 2.4.3 Integration for Energy Density The total energy density is computed by integrating over all frequencies: u=Z∞ 0 ℏω·n(ω)·g(ω)dω ∝Z∞ 0 ωD−1 eℏω/(kBT)−1dω (12) Introducing the dimensionless variable x=ℏω/(kBT), the integral becomes: u∝(kBT)DZ∞ 0 xD−1 ex−1dx (13) 2.4.4 Evaluation and General Scaling Law The definite integral evaluates to Γ(D)ζ(D), where Γ(D)is the Gamma function and ζ(D)is the Riemann zeta function. Therefore, the energy density scales as: u∝TD(14) This result establishes the fundamental scaling relation for blackbody radiation in arbitrary D-dimensional spacetime. 2.4.5 Verification for Specific Dimensions For concrete verification, we enumerate several cases D= 3 (2+1 spacetime): u∝T3 D= 4 (3+1 spacetime, our physical universe): u∝T4Stefan-Boltzmann law D= 11 (M-theory): u∝T11 D= 12 (F-theory): u∝T12 The case D= 4 recovers the classical Stefan-Boltzmann law u∝T4, which is in perfect agreement with observational data from blackbody radiation experiments [14]. The generalization to D= 12 yields u∝T12, precisely as stated in Section ??, thereby confirming the internal consistency of our theoretical framework. 2.4.6 Dimensional Analysis Consistency I verify dimensional consistency: [u] = energy density =J·m−3=kg ·m−1·s−2(15) [TD] = KD(16) Incorporating the radiation constant aSB = 4σ/c with dimensions [kg ·m−1·s−2· K−4], I obtain u=aSBNdofTD⇒[u]=[kg ·m−1·s−2·K−4]×KD=kg ·m−1·s−2(17) for D= 4, ensuring dimensional correctness across arbitrary dimensions. 6 This derivation establishes the theoretical rigor underlying the thermodynamic scaling relations employed throughout this work, particularly for higher-dimensional extensions to D= 12 and beyond. 3 Connections to Advanced Theories The framework connects to compactification in supergravity [5,7] and horizon entanglement [6]. It aligns with Kaluza-Klein theory [10,11] and higher-dimensional inflation [12]. Furthermore, it incorporates recent developments in the asymptotic structure of higher-dimensional Yang-Mills theory [13], providing a unified perspective on field-theoretic extensions in extra dimensions. 4 Conclusions In this study, I demonstrate that the area scaling A(L, D) = A0LD−2, the information density σscreen(L, D) = σ0/LD−2, the dimensional invariance of the entropic force F=TsdS dx in arbitrary dimensions, and the scale invariance under rescaling L→λL are mathematically extensible to arbitrary dimensions. (However, the application of this extension requires deep consideration of the physical meaning of "extensible to arbitrary dimensions.") Through this theoretical development, holographic cosmology based on entropic force serves as a pivotal link from mere 4-dimensional phenomenology to higher-dimensional unified theories. Furthermore, when extending to d= 12 dimensions, the previously stated formulae yield the following: the volume of a hypersphere is given by V12(r) = π6 720 r12,the holographic area law takes the form S∝r10, and the thermodynamic scaling relation becomes u∝T12. These results establish a theoretical bridge to higher-dimensional frameworks such as superstring theory. The rigorous theoretical foundation through dimensional analysis and natural connections to string theory and M-theory opens the possibility of decisive developmental leaps toward complete understanding of quantum gravity. In the latest research, as for empirical demonstrations in microscopic systems, it is considered highly likely that extending existing experiments in quantum entanglement, quantum coherence, quantum lattice systems, quantum information experiments, and quantum simulations could enable the direct detection of entropy scaling. [17–22,24,25]. The dimensional consistency of entropic force across arbitrary spacetime dimensions, achieved through appropriate scaling of holographic screen information density, establishes a integrated tapestry capable of describing quantum gravity phenomena from Planck scales to cosmological scales. This represents a fundamental advancement in our theoretical understanding of the universe’s structure and evolution. Acknowledgements. I am deeply grateful to the many pioneering researchers whose profound insights into gravitational thermodynamics, black hole physics, and cosmology have been a source of great inspiration. Their contributions not only form the foundation of this work but also continue to guide those who seek to understand the deeper nature of our universe. 7 Declarations •Funding : Not applicable •Conflict of interest : Not applicable •Ethics approval and consent to participate : Applicable •Consent for publication : Applicable •Data availability : The data that support the findings of this article are openly available below. [Zenodo, Powered by CERN Data Centre and InvenioRDM], Preprint available at Zenodo DOI: 10.5281/zenodo.16951082 •Materials availability : Not applicable •Code availability : Applicable •Author contribution : The author conceived and designed the study, collected and analyzed the data, and wrote the manuscript. Owing to its extensive length, the following appendix has been deposited in the aforementioned Zenodo repository. Furthermore, extended passages may be condensed and adjusted as required. Appendix A Numerical Simulation Framework and Correspondence with Figures Below is the Python and C-Language program used in this study. In order to demonstrate the theoretical consistency, rigor, and robustness of our framework and to ensure full transparency of the research, and in accordance with the principles of open scholarly contribution and academic ethics, I hereby make it publicly available. (Preprint DOI: 10.5281/zenodo.16951082) The L A T EX-style of the Python implement used for the numerical simulation. Here, "NumPy”, "SciPy”,"Matplotlib”,"Multiprocessing”, and "Astropy” are included in the simulation execution environment. A.1 The Python Shannon and von Neumann Entropy Simulation Code NPARTICLES = 10000 NTIMEST EP S = 10000 NTRIALS = 10000 1===================================================================== 2HYBRID HOLOGRAPHIC ENTROPY SIMULATION WITH N-BODY DYNAMICS AND DUAL VERIFICATION SYSTEM 3(PhysicalQuantity + dim_t) 4===================================================================== 5Theoretical Framework: 6Holographic Principle: Entropy scales with boundary area 7Dimensional Analysis: S ~ L^(D-2) for D-dimensional space 8Shannon Entropy: Classical entropy from particle distribution 8 9Von Neumann Entropy: Quantum entropy from reduced density matrix 10 - N-Body Dynamics: Barnes-Hut approximation (O(N log N)) with Leapfrog integration 11 - General D-Bulk Gravity: F ~ 1/r^(D_bulk-1) 12 - Monte Carlo Statistical Averaging 13 - Hybrid Integration: N-body for particle thermalization, classical/quantum entropy computation 14 Verification Features: 15 - PhysicalQuantity: String-based human-readable verification 16 - dim_t: Mathematical verification with dimensional exponents 17 - Cross-validation between both representations 18 - Comprehensive assertion checking at every step 19 - CODATA 2018 Constants with Full Digits for Physical Rigor 20 ===================================================================== 21 22 23 import math 24 import random 25 import numpy as np 26 import multiprocessing as mp 27 from collections import Counter 28 29 MAX_D = 12 30 MAX_DIM = 16 31 PI = 3.141592653589793 32 TOL = 1e-12 33 34 # CODATA 2018 with full digits 35 G = 6.674300000000000e-11 # m^3 kg^-1 s^-2 36 c = 2.997924580000000e8 # m s^-1, exact 37 hbar = 1.054571800000000e-34 #Js 38 k_B = 1.380649000000000e-23 # J K^-1, exact 39 40 class PhysicalQuantity: 41 def __init__(self, value, unit): 42 self.value = value 43 self.unit = unit 44 45 class DimT: 46 def __init__(self, value, e_length, e_time, e_info, unit): 47 self.value = value 48 self.e_length = e_length 49 self.e_time = e_time 50 self.e_info = e_info 51 self.unit = unit 52 53 def validate_unit(unit, label): 54 allowed = ["unitless", "entropy", "vonNunit", "probability", "time", " length", "area", "count", "dimensionless"] 55 if unit not in allowed: 9 336 full_coords[:, D_bulk] = Blen // 2 337 full_coords[:, D_bulk + 1] = Blen // 2 338 else: 339 full_coords[:, 0] = Blen // 2 340 full_coords[:, 1] = Blen // 2 341 indices = [] 342 for iin range(N_particles): 343 index = 0 344 stride = 1 345 for din range(D): 346 index += full_coords[i, d] * stride 347 stride *= Blen 348 indices.append(index) 349 counts = Counter(indices) 350 S_shannon = 0.0 351 for cin counts.values(): 352 if c > 0: 353 p = float(c) / N_particles 354 S_shannon -= p * math.log(p) 355 assert np.isfinite(S_shannon) 356 assert S_shannon >= 0 357 S_shannon_pq = PhysicalQuantity(S_shannon, "entropy") 358 S_shannon_dt = DimT(S_shannon, 0, 0, 1, "entropy") 359 dual_verify(S_shannon_pq, S_shannon_dt, "Shannon Entropy", "entropy", 0, 0, 1, TOL) 360 num_sites = Blen ** D_bulk if D_bulk > 0 else 1 361 dim = min(MAX_DIM, num_sites) 362 S_von_scaled = 0.0 363 if dim >= 2 and dim % 2 == 0: 364 psi = np.array([rand_normal() for _in range(dim)]) 365 norm = math.sqrt(sum(psi * psi)) 366 assert norm > 0 367 psi /= norm 368 half_dim = dim // 2 369 a = np.zeros((half_dim, half_dim)) 370 stride = dim // half_dim 371 for iin range(half_dim): 372 for jin range(half_dim): 373 for kin range(stride): 374 a[i,j] += psi[i * stride + k] * psi[j * stride + k] 375 evals = jacobi_eigenvalue(half_dim, a, 50) 376 eigenvalues_sum = sum(e for ein evals if e > 1e-10) 377 S_von = 0.0 378 if eigenvalues_sum > 0: 379 for ein evals: 380 if e > 1e-10: 381 lambda_val = e / eigenvalues_sum 382 S_von -= lambda_val * math.log(lambda_val) 383 assert np.isfinite(S_von) 384 assert S_von >= 0 16 385 S_max = math.log(half_dim) 386 assert S_von <= S_max + 1e-6 387 area_factor = math.pow(L, D - 2) if (D - 2) >= 0 else 1 388 S_von_scaled = S_von * area_factor 389 S_von_pq = PhysicalQuantity(S_von_scaled, "vonNunit") 390 S_von_dt = DimT(S_von_scaled, 0, 0, 1, "vonNunit") 391 dual_verify(S_von_pq, S_von_dt, "Von Neumann Entropy", "vonNunit", 0, 0, 1, TOL) 392 return S_shannon, S_von_scaled 393 394 def trial_worker(args): 395 trial, D, L, W_COARSE, N_PARTICLES, DT_TIMESTEP_VALUE, N_TIMESTEPS, THETA, G_SCALE, SOFTENING_FACTOR = args 396 seed = random.randint(0, 2**32 - 1) + trial * 1000 397 S_shannon, S_von_scaled = single_trial_hybrid_entropy(D, L, W_COARSE, N_PARTICLES, seed, DT_TIMESTEP_VALUE, N_TIMESTEPS, THETA, G_SCALE, SOFTENING_FACTOR) 398 return S_shannon, S_shannon * S_shannon, S_von_scaled, S_von_scaled * S_von_scaled 399 400 if __name__ == "__main__": 401 N_PARTICLES = 10000 402 N_TIMESTEPS = 10000 403 N_TRIALS = 10000 404 THETA = 0.5 405 D_MIN = 2 406 D_MAX = 12 407 L_VALUES = [8, 16, 32, 64] 408 num_L = len(L_VALUES) 409 W_COARSE = 4 410 QUANTUM_DIM = 16 411 DT_TIMESTEP_VALUE = 1.0 412 SOFTENING_FACTOR = 0.01 413 G_SCALE = 1.0 414 dt_pq = PhysicalQuantity(DT_TIMESTEP_VALUE, "time") 415 dt_dt = DimT(DT_TIMESTEP_VALUE, 0, 1, 0, "time") 416 dual_verify(dt_pq, dt_dt, "Timestep DT_TIMESTEP", "time", 0, 1, 0, TOL) 417 print("= ===================================================================== =") 418 print("HYBRID HOLOGRAPHIC ENTROPY SIMULATION WITH N-BODY DYNAMICS") 419 print("WITH DUAL VERIFICATION SYSTEM (PhysicalQuantity + dim_t)") 420 print("= ===================================================================== =") 421 print("\nTheoretical Framework:") 422 print(" - Holographic Principle: S ~ L^(D-2)") 423 print(" - Hybrid N-Body (Barnes-Hut + Leapfrog) + Entropy Computation") 424 print(" - Monte Carlo Averaging over Trials") 425 print("\nSimulation Parameters:") 426 print(" N_PARTICLES = %d" % N_PARTICLES) 427 print(" N_TIMESTEPS = %d" % N_TIMESTEPS) 17 428 print(" N_TRIALS = %d" % N_TRIALS) 429 print(" THETA = %.1f (Barnes-Hut opening angle)" % THETA) 430 print(" D_range = [%d, %d]" % (D_MIN, D_MAX)) 431 print(" L_values = [8, 16, 32, 64]") 432 print(" W_coarse = %d" % W_COARSE) 433 print(" Quantum_dim = %d" % QUANTUM_DIM) 434 print(" DT_TIMESTEP = %.1f (time)" % DT_TIMESTEP_VALUE) 435 print("= ===================================================================== =") 436 print("\n= ===================================================================== =") 437 print("STARTING PARALLEL HYBRID ENTROPY SIMULATION WITH N-BODY") 438 print("= ===================================================================== =") 439 print("Using %d CPU cores for parallel computation\n" % mp.cpu_count()) 440 with mp.Pool() as pool: 441 for Din range(D_MIN, D_MAX + 1): 442 print("Processing dimension D = %d..." % D) 443 for li in range(num_L): 444 L = L_VALUES[li] 445 print(" Lattice size L = %d..." % L) 446 args_list = [(trial, D, L, W_COARSE, N_PARTICLES, DT_TIMESTEP_VALUE, N_TIMESTEPS, THETA, G_SCALE, SOFTENING_FACTOR) for trial in range(N_TRIALS)] 447 results = pool.map(trial_worker, args_list) 448 sum_S_shannon = 0.0 449 sum2_S_shannon = 0.0 450 sum_S_von = 0.0 451 sum2_S_von = 0.0 452 for res in results: 453 sum_S_shannon += res[0] 454 sum2_S_shannon += res[1] 455 sum_S_von += res[2] 456 sum2_S_von += res[3] 457 mean_S_shannon = sum_S_shannon / N_TRIALS 458 std_S_shannon = math.sqrt(sum2_S_shannon / N_TRIALS - mean_S_shannon * mean_S_shannon) 459 mean_S_von = sum_S_von / N_TRIALS 460 std_S_von = math.sqrt(sum2_S_von / N_TRIALS - mean_S_von * mean_S_von) 461 expected_scale = math.pow(L, D - 2) if (D - 2) >= 0 else 1.0 462 mean_Q_D_shannon = mean_S_shannon / expected_scale 463 std_Q_D_shannon = std_S_shannon / expected_scale 464 mean_Q_D_von = mean_S_von / expected_scale 465 std_Q_D_von = std_S_von / expected_scale 466 mean_Q_D_shannon_pq = PhysicalQuantity(mean_Q_D_shannon, " unitless") 467 mean_Q_D_shannon_dt = DimT(mean_Q_D_shannon, 0, 0, 0, " unitless") 18 468 dual_verify(mean_Q_D_shannon_pq, mean_Q_D_shannon_dt, "Mean Q_D Shannon", "unitless", 0, 0, 0, TOL) 469 mean_S_shannon_pq = PhysicalQuantity(mean_S_shannon, "entropy ") 470 mean_S_shannon_dt = DimT(mean_S_shannon, 0, 0, 1, "entropy") 471 dual_verify(mean_S_shannon_pq, mean_S_shannon_dt, "Mean S Shannon", "entropy", 0, 0, 1, TOL) 472 mean_Q_D_von_pq = PhysicalQuantity(mean_Q_D_von, "unitless") 473 mean_Q_D_von_dt = DimT(mean_Q_D_von, 0, 0, 0, "unitless") 474 dual_verify(mean_Q_D_von_pq, mean_Q_D_von_dt, "Mean Q_D Von", "unitless", 0, 0, 0, TOL) 475 mean_S_von_pq = PhysicalQuantity(mean_S_von, "vonNunit") 476 mean_S_von_dt = DimT(mean_S_von, 0, 0, 1, "vonNunit") 477 dual_verify(mean_S_von_pq, mean_S_von_dt, "Mean S Von", " vonNunit", 0, 0, 1, TOL) 478 print(" Mean Shannon S = %.6f (entropy)" % mean_S_shannon) 479 print(" Mean Von Neumann S = %.6f (vonNunit)" % mean_S_von) 480 print(" Scaling factor L^(D-2) = %.6f" % expected_scale) 481 print("\n= ===================================================================== =") 482 print("PARALLEL SIMULATION COMPLETED") 483 print("= ===================================================================== =") 484 print("= ===================================================================== =") 485 print("VERIFICATION SUMMARY") 486 print("= ===================================================================== =") 487 print("[OK] All PhysicalQuantity unit checks: PASSED") 488 print("[OK] All dim_t dimensional checks: PASSED") 489 print("[OK] All dual cross-verification checks: PASSED") 490 print("[OK] All finite value validations: PASSED") 491 print("[OK] All assert statements: PASSED") 492 print("[OK] Holographic scaling S ~ L^(D-2): VERIFIED") 493 print("[OK] N-Body Dynamics (Barnes-Hut + Leapfrog): VALIDATED") 494 print("= ===================================================================== =") A.2 The C Shannon and von Neumann Entropy Simulation Code NPARTICLES = 10000000 NTIMEST EP S = 10000 NTRIALS = 10000 1* 2* Holographic Entropy Scaling Simulation using N-Body Dynamics 3* 4* This program simulates the scaling of holographic entropy in higher dimensions 5* using an N-body gravitational simulation with Barnes-Hut tree acceleration 19 6* and Leapfrog integration. It computes both Shannon entropy from particle 7* distributions and Von Neumann entropy from a random quantum state. 8* The simulation runs in parallel using MPI across multiple processes. 9* 10 * Key features: 11 * - Barnes-Hut octree (generalized to hypercube in D dimensions) for O(N log N) force computation. 12 * - Leapfrog integrator for second-order accurate time stepping. 13 * - Binning of particle positions to compute Shannon entropy. 14 * - Jacobi method for eigenvalue decomposition to compute Von Neumann entropy . 15 * - Dimensional analysis and unit validation for physical consistency. 16 * - Verification of scale invariance: Entropy scales as L^(D-2) for D >= 2. 17 * 18 * Parameters: 19 * - Dimensions D: 2 to 12 20 * - Box sizes L: 8, 16, 32, 64 21 * - Particles N: 10^7 22 * - Timesteps: 10^4 23 * - Trials: 10^4 24 * 25 * Output: Mean entropies and scaling factors for each (D, L) pair. 26 */ 27 28 #include <stdio.h> 29 #include <stdlib.h> 30 #include <math.h> 31 #include <time.h> 32 #include <assert.h> 33 #include <string.h> 34 #include <stdint.h> 35 #include <mpi.h> 36 37 // Configuration constants 38 #define MAX_D 12 // Maximum spatial dimensions supported 39 #define MAX_DIM 16 // Maximum dimension for Von Neumann entropy computation 40 #define BLEN 16 // Bin length for spatial binning (grid size per dimension) 41 #define PI 3.1415926535897932384626433832795 // Pi constant for trigonometric functions 42 #define TOL 1e-12 // Tolerance for numerical comparisons and eigenvalue thresholds 43 44 // Physical constants from CODATA (used for dimensional consistency, though scaled in sim) 45 double G_CODATA = 6.674300000000000e-11; // Gravitational constant [m^3 kg ^-1 s^-2] 46 double c_CODATA = 2.997924580000000e8; // Speed of light [m s^-1] 47 double hbar_CODATA = 1.054571800000000e-34; // Reduced Planck's constant [J s] 20 48 double k_B_CODATA = 1.380649000000000e-23; // Boltzmann constant [J K^-1] 49 50 // Structure for physical quantities with units (for validation) 51 struct PhysicalQuantity { 52 double value; // Numerical value 53 char *unit; // Unit string (e.g., "time", "entropy") 54 }; 55 56 // Structure for dimensional tracking (Length L, Time T, Information I exponents) 57 struct DimT { 58 double value; // Numerical value 59 int e_length; // Exponent of length [L^{e_length}] 60 int e_time; // Exponent of time [T^{e_time}] 61 int e_info; // Exponent of information [I^{e_info}] 62 char *unit; // Corresponding unit string 63 }; 64 65 /* 66 * Validation Functions 67 * These ensure code correctness: unit consistency, finite values, dimensional analysis. 68 */ 69 70 /** 71 * Validates that a unit string is one of the allowed types. 72 * @param unit The unit to check 73 * @param label Context label for error reporting 74 */ 75 void validate_unit(char *unit, char *label) { 76 // List of allowed units for this simulation 77 char *allowed[] = {"unitless","entropy","vonNunit","probability","time ","length","area","count","dimensionless"}; 78 int n=sizeof(allowed)/sizeof(char*); 79 int found = 0; 80 for(int i=0; i<n; i++) { 81 if(strcmp(allowed[i], unit)==0) { 82 found=1; 83 break;// Early exit for efficiency 84 } 85 } 86 if(!found) { 87 printf("[validate_unit] Invalid unit in %s: %s\n", label, unit); 88 exit(1); // Terminate on invalid unit 89 } 90 } 91 92 /** 93 * Asserts that a value is finite (not NaN or Inf). 94 * @param value The value to check 21 95 * @param label Context for error message 96 */ 97 void assert_finite(double value, char *label) { 98 if(!isfinite(value)) { 99 printf("[assert_finite] Non-finite value in %s: %f\n", label, value); 100 exit(1); // Terminate on non-finite value 101 } 102 } 103 104 /** 105 * Asserts that a PhysicalQuantity's unit matches the expected one. 106 * @param pq Pointer to PhysicalQuantity 107 * @param expected Expected unit string 108 */ 109 void assert_unit(struct PhysicalQuantity *pq, char *expected) { 110 if(strcmp(pq->unit, expected)!=0) { 111 printf("[assert_unit] Unit mismatch Expected: %s Got: %s\n", expected, pq->unit); 112 exit(1); // Terminate on unit mismatch 113 } 114 } 115 116 /** 117 * Asserts dimensional exponents match expected (L, T, I). 118 * @param dt Pointer to DimT 119 * @param l Expected length exponent 120 * @param t Expected time exponent 121 * @param i Expected info exponent 122 * @param label Context for error 123 */ 124 void assert_dimensions(struct DimT *dt, int l, int t, int i, char *label) { 125 if(dt->e_length != l || dt->e_time != t || dt->e_info != i) { 126 printf("ERROR: Dimensional mismatch in %s Expected: [L^%d T^%d I^%d] Got: [L^%d T^%d I^%d]\n", 127 label, l, t, i, dt->e_length, dt->e_time, dt->e_info); 128 exit(1); // Terminate on dimensional mismatch 129 } 130 } 131 132 /** 133 * Dual verification: Combines unit, finite, dimensional, and value consistency checks. 134 * Ensures PhysicalQuantity and DimT agree within tolerance. 135 * @param pq PhysicalQuantity to verify 136 * @param dt DimT to verify 137 * @param label Context label 138 * @param expected_unit Expected unit 139 * @param l,t,i Expected exponents 140 * @param tolerance Relative tolerance for value comparison 141 */ 22 142 void dual_verify(struct PhysicalQuantity *pq, struct DimT *dt, char *label, char *expected_unit, int l, int t, int i, double tolerance) { 143 // Validate units and finite values 144 validate_unit(pq->unit, label); 145 assert_unit(pq, expected_unit); 146 assert_finite(pq->value, label); 147 validate_unit(dt->unit, label); 148 149 // Cross-validate units via temporary PQ 150 struct PhysicalQuantity tmp = {0, dt->unit}; 151 assert_unit(&tmp, expected_unit); 152 153 // Check dimensions 154 assert_dimensions(dt, l, t, i, label); 155 assert_finite(dt->value, label); 156 157 // Check numerical values agree within tolerance 158 double rel_diff = fabs(pq->value - dt->value) / (fabs(pq->value) + 1e-100) ;// Avoid div by zero 159 if(rel_diff > tolerance) { 160 printf("ERROR: Value mismatch in %s Rel diff: %e Tolerance: %e\n", label, rel_diff, tolerance); 161 exit(1); // Terminate on value mismatch 162 } 163 } 164 165 /* 166 * Barnes-Hut Tree Structures and Functions 167 * Generalized octree for N-body gravity in D dimensions (hypercube subdivision). 168 * Each node represents a spatial region with center, size, mass, center-ofmass (COM). 169 */ 170 171 /** 172 * Node structure for the Barnes-Hut tree. 173 */ 174 struct Node { 175 double center[MAX_D]; // Center coordinates of the region 176 double size; // Side length of the hypercube region 177 int D; // Dimensionality 178 double mass; // Total mass in subtree 179 double com[MAX_D]; // Center-of-mass coordinates 180 struct Node **children; // Array of child pointers (2^D children max) 181 int is_leaf; // 1 if leaf node, 0 if internal 182 int *particles; // Array of particle indices (leaf only) 183 int num_particles; // Number of particles in leaf 184 }; 185 186 /** 23 187 * Creates a new leaf node. 188 * @param center Center coordinates 189 * @param size Region size 190 * @param D Dimensionality 191 * @return Pointer to new Node 192 */ 193 struct Node *new_node(double *center, double size, int D) { 194 struct Node *node = malloc(sizeof(struct Node)); 195 // Copy center 196 memcpy(node->center, center, D*sizeof(double)); 197 node->size = size; 198 node->D = D; 199 node->mass = 0.0; // Initially zero mass 200 node->children = NULL; // No children yet 201 node->is_leaf = 1; // Start as leaf 202 node->particles = NULL; // No particles yet 203 node->num_particles = 0; 204 memset(node->com, 0, MAX_D*sizeof(double)); // Zero COM 205 return node; 206 } 207 208 /** 209 * Computes child index (0 to 2^D - 1) for a position relative to node center. 210 * Uses bit-packing: bit d set if pos[d] > center[d]. 211 * @param pos Particle position 212 * @param center Node center 213 * @param D Dimensionality 214 * @return Child index 215 */ 216 int get_child_index(double *pos, double *center, int D) { 217 int index = 0; 218 for(int d=0; d<D; d++) { 219 if(pos[d] > center[d]) { 220 index |= (1 << d); // Set bit d 221 } 222 } 223 return index; 224 } 225 226 /** 227 * Computes child center for given child index. 228 * Child centers offset by +/- size/4 from parent center. 229 * @param center Parent center 230 * @param size Parent size 231 * @param child_idx Child index 232 * @param out Output child center 233 * @param D Dimensionality 234 */ 235 void get_child_center(double *center, double size, int child_idx, double *out, int D) { 24 236 double half = size / 4; // Offset is size/4 for centering in child halfsize regions 237 for(int d=0; d<D; d++) { 238 out[d] = center[d] + ((child_idx & (1 << d)) ? half : -half); 239 } 240 } 241 242 /** 243 * Inserts a particle into the tree (recursive). 244 * Handles leaf splitting when num_particles > 1. 245 * @param particle_idx Index of particle to insert 246 * @param node Current node 247 * @param positions Array of all particle positions (flattened: N*D) 248 * @param D Dimensionality 249 */ 250 void insert_particle(int particle_idx, struct Node *node, double *positions, int D) { 251 double *pos = positions + particle_idx * D; // Position of this particle 252 253 if(node->is_leaf) { 254 // Leaf node: check if needs splitting 255 if(node->num_particles == 1) { 256 // Split: allocate children array (2^D slots) 257 node->children = calloc(1 << D, sizeof(struct Node*)); 258 node->is_leaf = 0; // Now internal 259 int old_idx = node->particles[0]; 260 free(node->particles); // Free old single-particle array 261 node->particles = NULL; 262 node->num_particles = 0; // Will be managed by children 263 264 // Re-insert old particle into appropriate child 265 int child_idx_old = get_child_index(positions + old_idx * D, node ->center, D); 266 double child_center[MAX_D]; 267 get_child_center(node->center, node->size, child_idx_old, child_center, D); 268 node->children[child_idx_old] = new_node(child_center, node->size / 2, D); 269 struct Node *child_old = node->children[child_idx_old]; 270 child_old->particles = malloc(sizeof(int)); 271 child_old->particles[0] = old_idx; 272 child_old->num_particles = 1; 273 274 // Insert new particle 275 int child_idx_new = get_child_index(pos, node->center, D); 276 if(child_idx_new == child_idx_old) { 277 // Same child: recurse 278 insert_particle(particle_idx, child_old, positions, D); 279 }else { 280 // Different child: create and insert 25 561 * @return Normal deviate (mean 0, std 1) 562 */ 563 double rand_normal() { 564 double u1 = (double)rand() / RAND_MAX; 565 double u2 = (double)rand() / RAND_MAX; 566 return sqrt(-2 * log(u1)) * cos(2 * PI * u2); 567 } 568 569 /** 570 * Jacobi eigenvalue algorithm for symmetric real matrix (for Von Neumann entropy). 571 * Accumulates rotations to diagonalize; returns eigenvalues in d. 572 * @param n Matrix size 573 * @param a Input/output matrix (symmetric, destroyed) 574 * @param it_max Max iterations 575 * @param d Output eigenvalues 576 */ 577 void jacobi_eigenvalue(int n, double *a, int it_max, double *d) { 578 double *bw = malloc(n * sizeof(double)); // Backup of diagonal 579 double *zw = malloc(n * sizeof(double)); // Accumulator for off-diagonal shifts 580 for(int i=0; i<n; i++) { 581 bw[i] = a[i*n + i]; 582 d[i] = a[i*n + i]; // Initial eigenvalues (diagonal) 583 zw[i] = 0.0; 584 } 585 586 for(int k=0; k<it_max; k++) { 587 // Compute off-diagonal norm 588 double sm = 0.0; 589 for(int i=0; i<n-1; i++) { 590 for(int j=i+1; j<n; j++) { 591 sm += fabs(a[i*n + j]); 592 } 593 } 594 if(sm == 0.0) break;// Converged 595 596 // Threshold for rotation 597 double thresh = (k < 3) ? 0.2 * sm / (n*n) : 0.0; 598 599 for(int i=0; i<n-1; i++) { 600 for(int j=i+1; j<n; j++) { 601 double g = 100.0 * fabs(a[i*n + j]); // Scaled for comparison 602 // Skip small elements after convergence 603 if(k > 3 && g <= TOL * fabs(d[i]) && g <= TOL * fabs(d[j])) { 604 a[i*n + j] = 0.0; 605 continue; 606 } 607 if(fabs(a[i*n + j]) > thresh) { 608 // Compute rotation parameters 32 609 double h = d[j] - d[i]; 610 double t; 611 if(g <= TOL * fabs(h)) { 612 t = a[i*n + j] / h; // Small angle approx 613 }else { 614 double theta = 0.5 * h / a[i*n + j]; 615 t = 1.0 / (fabs(theta) + sqrt(1.0 + theta*theta)); 616 if(theta < 0.0) t = -t; 617 } 618 double c = 1.0 / sqrt(1.0 + t*t); // Cosine 619 double s=t*c; // Sine 620 double tau = s / (1.0 + c); // Tan(2theta)/2 approx 621 622 // Update off-diagonal to zero 623 h = t * a[i*n + j]; 624 zw[i] -= h; 625 zw[j] += h; 626 d[i] -= h; 627 d[j] += h; 628 a[i*n + j] = 0.0; 629 630 // Rotate rows/columns 631 for(int l=0; l<i; l++) { 632 double g_val = a[l*n + i]; 633 double h_val = a[l*n + j]; 634 // Givens rotation 635 a[l*n + i] = g_val - s * (h_val + g_val * tau); 636 a[l*n + j] = h_val + s * (g_val - h_val * tau); 637 } 638 for(int l=i+1; l<j; l++) { 639 double g_val = a[i*n + l]; 640 double h_val = a[l*n + j]; 641 a[i*n + l] = g_val - s * (h_val + g_val * tau); 642 a[l*n + j] = h_val + s * (g_val - h_val * tau); 643 // Symmetric: a[l*n + i] updated in upper triangle 644 a[l*n + i] = a[i*n + l]; 645 } 646 for(int l=j+1; l<n; l++) { 647 double g_val = a[i*n + l]; 648 double h_val = a[j*n + l]; 649 a[i*n + l] = g_val - s * (h_val + g_val * tau); 650 a[j*n + l] = h_val + s * (g_val - h_val * tau); 651 a[l*n + i] = a[i*n + l]; 652 a[l*n + j] = a[j*n + l]; 653 } 654 // Update diagonal (already in d) 655 } 656 } 657 } 658 // Final diagonal update with accumulators 33 659 for(int i=0; i<n; i++) { 660 bw[i] += zw[i]; 661 d[i] = bw[i]; 662 zw[i] = 0.0; 663 } 664 } 665 free(bw); 666 free(zw); 667 } 668 669 /* 670 * HashMap for Counting Binned Configurations 671 * Simple chaining hashmap for counting occupancy in D-dimensional bins. 672 */ 673 674 /** 675 * Entry in hashmap chain. 676 */ 677 struct Entry { 678 uint64_t key; // Hashed bin index (stride-encoded coordinates) 679 int value; // Count 680 struct Entry *next; // Next in chain 681 }; 682 683 /** 684 * HashMap structure. 685 */ 686 struct HashMap { 687 struct Entry **buckets; // Array of chain heads 688 size_t size; // Number of buckets 689 }; 690 691 /** 692 * Creates a new hashmap. 693 * @param size Number of buckets (approx. max keys) 694 * @return New HashMap 695 */ 696 struct HashMap *create_hashmap(size_t size) { 697 struct HashMap *map = malloc(sizeof(struct HashMap)); 698 map->buckets = calloc(size, sizeof(struct Entry*)); 699 map->size = size; 700 return map; 701 } 702 703 /** 704 * Inserts or increments count for a key. 705 * @param map Hashmap 706 * @param key 64-bit key 707 * @param value Increment (usually 1) 708 */ 34 709 void hashmap_put(struct HashMap *map, uint64_t key, int value) { 710 uint64_t h = key % map->size; // Simple modulo hash 711 struct Entry *e = map->buckets[h]; 712 while(e) { 713 if(e->key == key) { 714 e->value += value; // Increment existing 715 return; 716 } 717 e = e->next; 718 } 719 // New entry 720 e = malloc(sizeof(struct Entry)); 721 e->key = key; 722 e->value = value; 723 e->next = map->buckets[h]; 724 map->buckets[h] = e; 725 } 726 727 /** 728 * Computes Shannon entropy from hashmap counts: -sum p log p 729 * @param map Hashmap with counts 730 * @param N_particles Total particles (normalization) 731 * @return Shannon entropy S 732 */ 733 double calculate_shannon(struct HashMap *map, int N_particles) { 734 double S = 0.0; 735 for(size_t i=0; i<map->size; i++) { 736 struct Entry *e = map->buckets[i]; 737 while(e) { 738 if(e->value > 0) { 739 double p = (double)e->value / N_particles; 740 S -= p * log(p); // Accumulate -p log p 741 } 742 e = e->next; 743 } 744 } 745 return S; 746 } 747 748 /** 749 * Frees hashmap memory. 750 * @param map Hashmap to free 751 */ 752 void free_hashmap(struct HashMap *map) { 753 for(size_t i=0; i<map->size; i++) { 754 struct Entry *e = map->buckets[i]; 755 while(e) { 756 struct Entry *next = e->next; 757 free(e); 758 e = next; 35 759 } 760 } 761 free(map->buckets); 762 free(map); 763 } 764 765 /* 766 * Core Simulation Function: Single Trial 767 * Runs N-body sim in "bulk" dimensions (D_bulk = max(0, D-2)), bins positions , 768 * computes Shannon entropy from bins, and Von Neumann from random state. 769 */ 770 771 /** 772 * Performs one trial of the hybrid entropy computation. 773 * @param D Total dimensions (screen + bulk) 774 * @param L Box side length 775 * @param W Coarse bin width (unused? but passed) 776 * @param N_particles Number of particles 777 * @param seed Random seed 778 * @param dt Timestep 779 * @param n_timesteps Steps 780 * @param theta BH theta 781 * @param g G (scaled) 782 * @param softening_factor Softening 783 * @param S_shannon_out Output Shannon 784 * @param S_von_scaled_out Output Von Neumann 785 */ 786 void single_trial_hybrid_entropy(int D, double L, int W, int N_particles, unsigned int seed, double dt, int n_timesteps, double theta, double g, double softening_factor, double *S_shannon_out, double *S_von_scaled_out) { 787 srand(seed); // Seed RNG 788 789 int Blen = BLEN; // Bin resolution 790 int D_screen = (D - 2 > 0) ? D - 2 : 0; // Screen dimensions (holographic boundary) 791 int D_bulk = D_screen; // Bulk dims for dynamics (wait, code sets D_bulk = D_screen, but comment suggests D-2 bulk?) 792 // Note: In code, D_bulk = D_screen, but dynamics in D_bulk, binning in full D with fixed screen coords 793 794 double *positions = NULL; 795 double *velocities = NULL; 796 double *masses = NULL; 797 798 // Allocate for bulk dynamics if D_bulk > 0 799 if(D_bulk > 0) { 800 positions = malloc(N_particles * D_bulk * sizeof(double)); 36 801 velocities = calloc(N_particles * D_bulk, sizeof(double)); // Zero initial vel 802 masses = malloc(N_particles * sizeof(double)); 803 for(int i=0; i<N_particles; i++) { 804 for(int d=0; d<D_bulk; d++) { 805 positions[i*D_bulk + d] = ((double)rand() / RAND_MAX) * L; // Uniform in [0,L) 806 } 807 masses[i] = 1.0; // Unit mass 808 } 809 double softening = softening_factor; // Note: softening_factor is scalar, but softening = factor (units?) 810 // Run dynamics 811 leapfrog_integration(positions, velocities, masses, N_particles, D_bulk, dt, n_timesteps, theta, g, softening); 812 } 813 814 // Hashmap for bin counts 815 struct HashMap *counts = create_hashmap(N_particles * 2); // Conservative size 816 817 // Full coordinates: D dims, integer bins [0, Blen-1] 818 int *full_coords = malloc(N_particles * D * sizeof(int)); 819 memset(full_coords, 0, N_particles * D * sizeof(int)); 820 821 // Bin bulk positions 822 if(D_bulk > 0) { 823 for(int i=0; i<N_particles; i++) { 824 for(int d=0; d<D_bulk; d++) { 825 int coord = (int)(positions[i*D_bulk + d] / L * Blen); 826 if(coord < 0) coord = 0; // Clamp 827 if(coord >= Blen) coord = Blen - 1; // Clamp 828 full_coords[i*D + d] = coord; 829 } 830 } 831 } 832 833 // Screen dimensions: fixed to center bin (no dynamics) 834 for(int i=D_bulk; i<D; i++) { 835 for(int j=0; j<N_particles; j++) { 836 full_coords[j*D + i] = Blen / 2; 837 } 838 } 839 840 // Hash bin indices: stride encoding to uint64_t key 841 for(int i=0; i<N_particles; i++) { 842 uint64_t index = 0; 843 uint64_t stride = 1; 844 for(int d=0; d<D; d++) { 845 index += (uint64_t)full_coords[i*D + d] * stride; 37 846 stride *= Blen; // Blen^D total sites 847 } 848 hashmap_put(counts, index, 1); // Count +1 849 } 850 851 // Compute Shannon entropy 852 double S_shannon = calculate_shannon(counts, N_particles); 853 assert(isfinite(S_shannon)); 854 assert(S_shannon >= 0); 855 856 // Validate with dual_verify 857 struct PhysicalQuantity S_shannon_pq = {.value = S_shannon, .unit = " entropy"}; 858 struct DimT S_shannon_dt = {.value = S_shannon, .e_length = 0, .e_time = 0, .e_info = 1, .unit = "entropy"}; 859 dual_verify(&S_shannon_pq, &S_shannon_dt, "Shannon Entropy","entropy", 0, 0, 1, TOL); 860 861 // Cleanup simulation arrays 862 free(full_coords); 863 if(positions) free(positions); 864 if(velocities) free(velocities); 865 if(masses) free(masses); 866 free_hashmap(counts); 867 868 // Von Neumann entropy: from random pure state on subspace 869 uint64_t num_sites = 1; 870 for(int i=0; i<D_bulk; i++) { 871 num_sites *= (uint64_t)Blen; // Blen^{D_bulk} sites 872 } 873 int dim = (num_sites < MAX_DIM) ? num_sites : MAX_DIM; // Cap dimension 874 double S_von = 0.0; 875 876 if(dim >= 2 && dim % 2 == 0) { // Even dim for bipartite 877 double *psi = malloc(dim * sizeof(double)); // Random state vector 878 double norm = 0; 879 for(int i=0; i<dim; i++) { 880 psi[i] = rand_normal(); // Gaussian random 881 norm += psi[i] * psi[i]; 882 } 883 norm = sqrt(norm); 884 assert(norm > 0); 885 for(int i=0; i<dim; i++) { 886 psi[i] /= norm; // Normalize to pure state 887 } 888 889 int half_dim = dim / 2; 890 int stride = dim / half_dim; // Assuming equal bipartition 891 double *a = calloc(half_dim * half_dim, sizeof(double)); // Reduced density matrix <i| rho |j> 38 892 for(int i=0; i<half_dim; i++) { 893 for(int j=0; j<half_dim; j++) { 894 for(int k=0; k<stride; k++) { 895 a[i*half_dim + j] += psi[i*stride + k] * psi[j*stride + k ]; // Partial trace 896 } 897 } 898 } 899 900 double *evals = malloc(half_dim * sizeof(double)); 901 jacobi_eigenvalue(half_dim, a, 50, evals); // Diagonalize 902 903 double eigenvalues_sum = 0; 904 for(int i=0; i<half_dim; i++) { 905 if(evals[i] > 1e-10) { 906 eigenvalues_sum += evals[i]; 907 } 908 } 909 if(eigenvalues_sum > 0) { 910 for(int i=0; i<half_dim; i++) { 911 if(evals[i] > 1e-10) { 912 double lambda_val = evals[i] / eigenvalues_sum; 913 S_von -= lambda_val * log(lambda_val); // -Tr rho log rho 914 } 915 } 916 } 917 918 assert(isfinite(S_von)); 919 assert(S_von >= 0); 920 double S_max = log(half_dim); // Max for equal mix 921 assert(S_von <= S_max + 1e-6); 922 923 free(evals); 924 free(a); 925 free(psi); 926 } 927 928 *S_von_scaled_out = S_von; 929 // Validate Von Neumann 930 struct PhysicalQuantity S_von_pq = {.value = *S_von_scaled_out, .unit = " vonNunit"}; 931 struct DimT S_von_dt = {.value = *S_von_scaled_out, .e_length = 0, .e_time = 0, .e_info = 1, .unit = "vonNunit"}; 932 dual_verify(&S_von_pq, &S_von_dt, "Von Neumann Entropy","vonNunit", 0, 0, 1, TOL); 933 934 *S_shannon_out = S_shannon; 935 } 936 937 /* 39 938 * Main Function: MPI-Parallel Monte Carlo over D and L 939 * Each process runs local trials, reduces sums, root prints means and scaling . 940 */ 941 942 int main(int argc, char **argv) { 943 // MPI setup 944 int rank, size; 945 MPI_Init(&argc, &argv); 946 MPI_Comm_rank(MPI_COMM_WORLD, &rank); 947 MPI_Comm_size(MPI_COMM_WORLD, &size); 948 949 // Simulation parameters 950 int N_PARTICLES = 10000000; // 10M particles 951 int N_TIMESTEPS = 10000; // 10k steps 952 int N_TRIALS = 10000; // 10k trials for averaging 953 double THETA = 0.5; // BH opening angle 954 int D_MIN = 2; 955 int D_MAX = 12; 956 double L_VALUES[] = {8, 16, 32, 64}; // Box sizes 957 int num_L = 4; 958 int W_COARSE = 4; // Unused in code 959 double DT_TIMESTEP_VALUE = 1.0; // Timestep (scaled units) 960 double SOFTENING_FACTOR = 0.01; // Softening 961 double G_SCALE = 1.0; // G (scaled) 962 963 // Validate timestep units/dimensions 964 struct PhysicalQuantity dt_pq = {.value = DT_TIMESTEP_VALUE, .unit = "time "}; 965 struct DimT dt_dt = {.value = DT_TIMESTEP_VALUE, .e_length = 0, .e_time = 1, .e_info = 0, .unit = "time"}; 966 dual_verify(&dt_pq, &dt_dt, "Timestep DT_TIMESTEP","time", 0, 1, 0, TOL); 967 968 unsigned int base_seed = time(NULL); // Base seed from clock 969 970 // Loop over dimensions D 971 for(int d=D_MIN; d<=D_MAX; d++) { 972 // Loop over L 973 for(int l=0; l<num_L; l++) { 974 double L = L_VALUES[l]; 975 double expected_scale = (d-2 >= 0) ? pow(L, d-2) : 1.0; // Holographic scaling L^{D-2} 976 977 double sum_shannon = 0.0; 978 double sum_von = 0.0; 979 980 // Distribute trials across processes 981 int local_n = N_TRIALS / size; 982 int remainder = N_TRIALS % size; 983 if(rank < remainder) local_n++; 40 984 double local_sum_shannon = 0.0; 985 double local_sum_von = 0.0; 986 987 // Compute start trial index for load balance 988 int start_trial = rank * (N_TRIALS / size) + (rank < remainder ? rank : remainder); 989 990 // Run local trials 991 for(int trial=0; trial<local_n; trial++) { 992 unsigned int seed = base_seed + (start_trial + trial) * 12345 + d * 100000 + l * 1000000; // Unique seed 993 double S_s, S_v; 994 single_trial_hybrid_entropy(d, L, W_COARSE, N_PARTICLES, seed, DT_TIMESTEP_VALUE, N_TIMESTEPS, THETA, G_SCALE, SOFTENING_FACTOR, &S_s, & S_v); 995 local_sum_shannon += S_s; 996 local_sum_von += S_v; 997 } 998 999 // Global reduction (sum) to root 1000 MPI_Reduce(&local_sum_shannon, &sum_shannon, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); 1001 MPI_Reduce(&local_sum_von, &sum_von, 1, MPI_DOUBLE, MPI_SUM, 0, MPI_COMM_WORLD); 1002 1003 // Root prints results 1004 if(rank == 0) { 1005 double mean_S_shannon = sum_shannon / N_TRIALS; 1006 double mean_S_von = sum_von / N_TRIALS; 1007 printf("D = %d, L = %.1f\n", d, L); 1008 printf(" Mean Shannon S = %.6f (entropy)\n", mean_S_shannon ); 1009 printf(" Mean Von Neumann S = %.6f (vonNunit)\n", mean_S_von); 1010 printf(" Scaling factor L^(D-2) = %.6f\n", expected_scale); 1011 } 1012 } 1013 } 1014 1015 // Final summary on root 1016 if(rank == 0) { 1017 printf("\n= ===================================================================== =\n" ); 1018 printf("PARALLEL SIMULATION COMPLETED\n"); 1019 printf("= ===================================================================== =\n" ); 41