1#[allow(unused_imports)]
6use super::functions::*;
7
8#[cfg(test)]
9mod tests {
10 use super::*;
11 use crate::fluid_sim_gpu::FlipParticle;
12 use crate::fluid_sim_gpu::FluidSimConfig;
13 use crate::fluid_sim_gpu::FluidSimulation;
14 use crate::fluid_sim_gpu::GpuBoundaryBox;
15 use crate::fluid_sim_gpu::GpuNeighborList;
16 use crate::fluid_sim_gpu::LbmCellType;
17 use crate::fluid_sim_gpu::LbmD2Q9;
18 use crate::fluid_sim_gpu::MacGrid;
19 use crate::fluid_sim_gpu::MultiGpuDomain;
20 use crate::fluid_sim_gpu::SphConfig;
21 use crate::fluid_sim_gpu::SphKernels;
22 use crate::fluid_sim_gpu::SphParticle;
23 #[test]
24 fn test_poly6_inside() {
25 let w = SphKernels::poly6(0.0, 0.1);
26 assert!(w > 0.0);
27 }
28 #[test]
29 fn test_poly6_outside() {
30 let w = SphKernels::poly6(0.2, 0.1);
31 assert!((w).abs() < 1e-15);
32 }
33 #[test]
34 fn test_spiky_grad_inside() {
35 let r_vec = [0.05, 0.0, 0.0];
36 let g = SphKernels::spiky_grad(r_vec, 0.05, 0.1);
37 assert!(g[0] != 0.0);
38 }
39 #[test]
40 fn test_viscosity_laplacian() {
41 let lap = SphKernels::viscosity_laplacian(0.05, 0.1);
42 assert!(lap > 0.0);
43 }
44 #[test]
45 fn test_sph_compute_density_single_particle() {
46 let mut particles = vec![SphParticle::new([0.0; 3], 1.0)];
47 let config = SphConfig::default();
48 sph_compute_density(&mut particles, &config);
49 assert!(particles[0].density > 0.0);
50 }
51 #[test]
52 fn test_sph_step_does_not_panic() {
53 let mut particles = vec![
54 SphParticle::new([0.0, 0.0, 0.0], 0.001),
55 SphParticle::new([0.05, 0.0, 0.0], 0.001),
56 ];
57 let config = SphConfig::default();
58 sph_step(&mut particles, &config);
59 }
60 #[test]
61 fn test_lbm_d2q9_new() {
62 let lbm = LbmD2Q9::new(8, 8, 1.0);
63 assert_eq!(lbm.nx, 8);
64 assert_eq!(lbm.ny, 8);
65 }
66 #[test]
67 fn test_lbm_density_initial() {
68 let lbm = LbmD2Q9::new(4, 4, 1.0);
69 let rho = lbm.density(2, 2);
70 assert!((rho - 1.0).abs() < 1e-10);
71 }
72 #[test]
73 fn test_lbm_velocity_initial_zero() {
74 let lbm = LbmD2Q9::new(4, 4, 1.0);
75 let [ux, uy] = lbm.velocity(2, 2);
76 assert!(ux.abs() < 1e-10);
77 assert!(uy.abs() < 1e-10);
78 }
79 #[test]
80 fn test_lbm_f_equilibrium_rest() {
81 let f_eq = LbmD2Q9::f_equilibrium(1.0, 0.0, 0.0, 0);
82 assert!((f_eq - D2Q9_W[0]).abs() < 1e-12);
83 }
84 #[test]
85 fn test_lbm_step_does_not_panic() {
86 let mut lbm = LbmD2Q9::new(8, 8, 0.6);
87 lbm.set_inlet_velocity(0.05, 0.0);
88 lbm.step();
89 }
90 #[test]
91 fn test_lbm_bounce_back_solid() {
92 let mut lbm = LbmD2Q9::new(4, 4, 1.0);
93 lbm.cell_type[0] = LbmCellType::Solid;
94 lbm.apply_bounce_back();
95 }
96 #[test]
97 fn test_mac_grid_new() {
98 let g = MacGrid::new(4, 4, 4, 0.1);
99 assert_eq!(g.nx, 4);
100 assert_eq!(g.u.len(), 5 * 4 * 4);
101 }
102 #[test]
103 fn test_mac_grid_set_get_u() {
104 let mut g = MacGrid::new(4, 4, 4, 0.1);
105 g.set_u(2, 1, 1, 3.125);
106 assert!((g.get_u(2, 1, 1) - 3.125).abs() < 1e-12);
107 }
108 #[test]
109 fn test_mac_grid_compute_divergence_zero() {
110 let mut g = MacGrid::new(4, 4, 4, 0.1);
111 g.flags = vec![1u8; 4 * 4 * 4];
112 g.compute_divergence();
113 for &d in &g.div {
114 assert!(d.abs() < 1e-12);
115 }
116 }
117 #[test]
118 fn test_jacobi_pressure_solve() {
119 let mut g = MacGrid::new(4, 4, 4, 0.1);
120 g.flags = vec![1u8; 4 * 4 * 4];
121 g.jacobi_pressure_solve(1000.0, 0.01, 5);
122 }
123 #[test]
124 fn test_compute_vorticity_zero_vel() {
125 let g = MacGrid::new(8, 8, 8, 0.1);
126 let vort = compute_vorticity(&g);
127 assert_eq!(vort.len(), 8 * 8 * 8);
128 for v in &vort {
129 assert!(v[0].abs() < 1e-12);
130 }
131 }
132 #[test]
133 fn test_surface_tension_csf() {
134 let nx = 4;
135 let ny = 4;
136 let nz = 4;
137 let phi = vec![-1.0f64; nx * ny * nz];
138 let forces = surface_tension_csf(&phi, nx, ny, nz, 0.1, 0.0728);
139 assert_eq!(forces.len(), nx * ny * nz);
140 }
141 #[test]
142 fn test_flip_particle_new() {
143 let p = FlipParticle::new([1.0, 2.0, 3.0]);
144 assert_eq!(p.position, [1.0, 2.0, 3.0]);
145 assert_eq!(p.velocity, [0.0; 3]);
146 }
147 #[test]
148 fn test_p2g_g2p_roundtrip() {
149 let mut particles = vec![FlipParticle::new([0.25, 0.25, 0.25])];
150 particles[0].velocity = [1.0, 0.0, 0.0];
151 let mut grid = MacGrid::new(4, 4, 4, 0.1);
152 p2g_transfer(&particles, &mut grid);
153 }
154 #[test]
155 fn test_fluid_simulation_new() {
156 let config = FluidSimConfig::default();
157 let sim = FluidSimulation::new(config);
158 assert_eq!(sim.particles.len(), 0);
159 assert!((sim.time).abs() < 1e-15);
160 }
161 #[test]
162 fn test_fluid_simulation_add_block() {
163 let config = FluidSimConfig::default();
164 let mut sim = FluidSimulation::new(config);
165 sim.add_fluid_block([0.1; 3], [0.9; 3], 8);
166 assert!(sim.particles.len() >= 8);
167 }
168 #[test]
169 fn test_fluid_simulation_step() {
170 let config = FluidSimConfig {
171 grid_size: [8, 8, 8],
172 pressure_iters: 5,
173 ..Default::default()
174 };
175 let mut sim = FluidSimulation::new(config);
176 sim.add_fluid_block([0.1; 3], [0.5; 3], 4);
177 sim.step();
178 assert_eq!(sim.step_count, 1);
179 assert!(sim.time > 0.0);
180 }
181 #[test]
182 fn test_fluid_kinetic_energy_increases_under_gravity() {
183 let config = FluidSimConfig {
184 grid_size: [8, 8, 8],
185 pressure_iters: 5,
186 ..Default::default()
187 };
188 let mut sim = FluidSimulation::new(config);
189 sim.add_fluid_block([0.4; 3], [0.6; 3], 4);
190 let ke0 = sim.kinetic_energy();
191 sim.step();
192 let ke1 = sim.kinetic_energy();
193 assert!(ke1 >= ke0);
194 }
195 #[test]
196 fn test_advect_scalar_uniform_field() {
197 let nx = 4;
198 let ny = 4;
199 let nz = 4;
200 let grid = MacGrid::new(nx, ny, nz, 0.1);
201 let phi = vec![1.0f64; nx * ny * nz];
202 let result = advect_scalar(&phi, &grid, 0.01, [0.0, -9.81, 0.0]);
203 assert_eq!(result.len(), nx * ny * nz);
204 for v in &result {
205 assert!((v - 1.0).abs() < 1e-10);
206 }
207 }
208 #[test]
209 fn test_sph_particle_gravity_descent() {
210 let mut p = vec![SphParticle::new([0.0, 1.0, 0.0], 0.001)];
211 let config = SphConfig {
212 gravity: [0.0, -9.81, 0.0],
213 ..Default::default()
214 };
215 let y0 = p[0].position[1];
216 for _ in 0..100 {
217 sph_step(&mut p, &config);
218 }
219 assert!(p[0].position[1] < y0);
220 }
221 #[test]
222 fn test_lbm_outlet_pressure() {
223 let mut lbm = LbmD2Q9::new(8, 8, 0.8);
224 lbm.set_outlet_pressure(1.0);
225 assert_eq!(lbm.cell_type[7], LbmCellType::Outlet);
226 }
227 #[test]
228 fn test_mac_grid_pressure_project_no_panic() {
229 let mut g = MacGrid::new(4, 4, 4, 0.1);
230 g.flags = vec![1u8; 64];
231 g.pressure_project(1000.0, 0.01);
232 }
233 #[test]
234 fn test_gpu_sph_density_summation_parallel() {
235 let mut particles = vec![
236 SphParticle::new([0.0, 0.0, 0.0], 0.01),
237 SphParticle::new([0.05, 0.0, 0.0], 0.01),
238 ];
239 let config = SphConfig::default();
240 gpu_sph_density_parallel(&mut particles, &config);
241 assert!(particles[0].density > 0.0);
242 assert!(particles[1].density > 0.0);
243 }
244 #[test]
245 fn test_gpu_pressure_jacobi_single_iter() {
246 let mut g = MacGrid::new(4, 4, 4, 0.1);
247 g.flags = vec![1u8; 64];
248 gpu_jacobi_pressure_solve(&mut g, 1000.0, 0.01, 3);
249 assert_eq!(g.p.len(), 64);
250 }
251 #[test]
252 fn test_gpu_neighbor_list_update() {
253 let positions = vec![[0.0f64, 0.0, 0.0], [0.05, 0.0, 0.0], [0.3, 0.0, 0.0]];
254 let nl = GpuNeighborList::build(&positions, 0.1, [0.5, 0.5, 0.5]);
255 let neighbors = nl.neighbors_of(0);
256 assert!(neighbors.contains(&1));
257 assert!(!neighbors.contains(&2));
258 }
259 #[test]
260 fn test_gpu_lbm_bgk_collision() {
261 let mut lbm = LbmD2Q9::new(8, 8, 0.8);
262 gpu_lbm_bgk_collide(&mut lbm);
263 let rho = lbm.density(4, 4);
264 assert!((rho - 1.0).abs() < 1e-9);
265 }
266 #[test]
267 fn test_morton_sort_particles() {
268 let mut particles = vec![
269 SphParticle::new([0.9, 0.9, 0.9], 1.0),
270 SphParticle::new([0.0, 0.0, 0.0], 1.0),
271 SphParticle::new([0.5, 0.5, 0.5], 1.0),
272 ];
273 morton_sort_particles(&mut particles, [1.0, 1.0, 1.0]);
274 assert!(particles[0].position[0] < 0.2);
275 }
276 #[test]
277 fn test_gpu_particle_integrate_euler() {
278 let mut particles = vec![SphParticle::new([0.0, 0.0, 0.0], 1.0)];
279 particles[0].force = [1.0, 0.0, 0.0];
280 particles[0].density = 1.0;
281 gpu_particle_integrate_euler(&mut particles, 0.01);
282 assert!(particles[0].velocity[0] > 0.0);
283 assert!(particles[0].position[0] > 0.0);
284 }
285 #[test]
286 fn test_gpu_particle_integrate_verlet() {
287 let mut particles = vec![SphParticle::new([0.0, 1.0, 0.0], 1.0)];
288 particles[0].force = [0.0, -9.81, 0.0];
289 particles[0].density = 1.0;
290 let dt = 0.01;
291 let y0 = particles[0].position[1];
292 gpu_particle_integrate_verlet(&mut particles, dt, dt);
293 assert!(particles[0].position[1] < y0);
294 }
295 #[test]
296 fn test_gpu_boundary_condition_box() {
297 let mut particles = vec![SphParticle::new([1.5, 0.5, 0.5], 1.0)];
298 particles[0].velocity = [1.0, 0.0, 0.0];
299 let bounds = GpuBoundaryBox {
300 min: [0.0; 3],
301 max: [1.0; 3],
302 restitution: 0.5,
303 };
304 gpu_apply_boundary_box(&mut particles, &bounds);
305 assert!(particles[0].position[0] <= 1.0);
306 assert!(particles[0].velocity[0] <= 0.0);
307 }
308 #[test]
309 fn test_multi_gpu_domain_decomp_2_devices() {
310 let positions: Vec<[f64; 3]> = (0..8).map(|i| [i as f64 * 0.1, 0.0, 0.0]).collect();
311 let domains = MultiGpuDomain::decompose_x(&positions, 2);
312 assert_eq!(domains.len(), 2);
313 let total: usize = domains.iter().map(|d| d.particle_indices.len()).sum();
314 assert_eq!(total, 8);
315 }
316 #[test]
317 fn test_gpu_reduce_energy() {
318 let particles = [
319 SphParticle::new([0.0; 3], 1.0),
320 SphParticle::new([0.0; 3], 1.0),
321 ];
322 let mut p0 = particles[0].clone();
323 let mut p1 = particles[1].clone();
324 p0.velocity = [2.0, 0.0, 0.0];
325 p1.velocity = [0.0, 1.0, 0.0];
326 let particles2 = vec![p0, p1];
327 let ke = gpu_reduce_kinetic_energy(&particles2, 0.5);
328 assert!((ke - 1.25).abs() < 1e-10);
329 }
330 #[test]
331 fn test_gpu_reduce_momentum() {
332 let mut p = SphParticle::new([0.0; 3], 1.0);
333 p.velocity = [3.0, 0.0, 0.0];
334 let particles = vec![p];
335 let mom = gpu_reduce_momentum(&particles, 2.0);
336 assert!((mom[0] - 6.0).abs() < 1e-10);
337 }
338 #[test]
339 fn test_gpu_advect_field() {
340 let nx = 4usize;
341 let ny = 4;
342 let field: Vec<f64> = (0..nx * ny).map(|i| i as f64).collect();
343 let vel = vec![[0.0f64; 2]; nx * ny];
344 let result = gpu_advect_2d(&field, &vel, nx, ny, 0.1, 0.01);
345 assert_eq!(result.len(), nx * ny);
346 }
347 #[test]
348 fn test_gpu_pressure_poisson_staggered() {
349 let mut p_grid = vec![0.0f64; 4 * 4];
350 let div = vec![0.0f64; 4 * 4];
351 gpu_pressure_poisson_jacobi_2d(&mut p_grid, &div, 4, 4, 0.1, 5);
352 for &v in &p_grid {
353 assert!(v.abs() < 1e-10);
354 }
355 }
356 #[test]
357 fn test_neighbor_list_cell_counts() {
358 let positions: Vec<[f64; 3]> = (0..4).map(|i| [i as f64 * 0.02, 0.0, 0.0]).collect();
359 let nl = GpuNeighborList::build(&positions, 0.1, [0.5, 0.5, 0.5]);
360 let n0 = nl.neighbors_of(0);
361 assert_eq!(n0.len(), 3);
362 }
363 #[test]
364 fn test_morton_code_ordering() {
365 let a = morton_encode_3d(0, 0, 0);
366 let b = morton_encode_3d(1, 0, 0);
367 let c = morton_encode_3d(0, 1, 0);
368 assert!(a < b);
369 assert!(a < c);
370 }
371 #[test]
372 fn test_gpu_sph_verlet_energy_conservation_approx() {
373 let mut particles = vec![SphParticle::new([0.0, 0.5, 0.0], 0.01)];
374 let config = SphConfig {
375 gravity: [0.0, 0.0, 0.0],
376 ..Default::default()
377 };
378 sph_compute_density(&mut particles, &config);
379 sph_compute_pressure(&mut particles, &config);
380 let v_before = particles[0].velocity;
381 gpu_particle_integrate_verlet(&mut particles, 0.001, 0.001);
382 let v_after = particles[0].velocity;
383 assert!((v_before[1] - v_after[1]).abs() < 1e-8);
384 }
385 #[test]
386 fn test_domain_decomp_particles_go_to_correct_half() {
387 let positions: Vec<[f64; 3]> = vec![
388 [0.1, 0.0, 0.0],
389 [0.2, 0.0, 0.0],
390 [0.6, 0.0, 0.0],
391 [0.8, 0.0, 0.0],
392 ];
393 let domains = MultiGpuDomain::decompose_x(&positions, 2);
394 assert_eq!(domains[0].particle_indices.len(), 2);
395 assert_eq!(domains[1].particle_indices.len(), 2);
396 }
397}