PySplashsurf: f32 and f64 variants for tests

pysplashsurf/tests/test_basic.py

@@ -108,8 +108,8 @@ def test_pipeline_f64():
     impl_basic_test(np.float64)
 
 
-def test_reconstruct():
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+def reconstruct_test(dtype):
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
 
     reconstruction = pysplashsurf.reconstruct_surface(
         particles,
@@ -130,17 +130,25 @@ def test_reconstruct():
 
     mesh = reconstruction.mesh
 
-    assert mesh.dtype == np.float32
+    assert mesh.dtype == dtype
 
-    assert reconstruction.particle_densities.dtype == np.float32
+    assert reconstruction.particle_densities.dtype == dtype
     assert len(reconstruction.particle_densities) == len(particles)
 
     assert len(mesh.vertices) in range(25000, 30000)
     assert len(mesh.triangles) in range(49000, 53000)
 
 
-def test_neighborhood_search():
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+def test_reconstruct_f32():
+    reconstruct_test(np.float32)
+
+
+def test_reconstruct_f64():
+    reconstruct_test(np.float64)
+
+
+def neighborhood_search_test(dtype):
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
 
     reconstruction = pysplashsurf.reconstruct_surface(
         particles,
@@ -174,8 +182,16 @@ def test_neighborhood_search():
     # TODO: Compare with naive neighbor search
 
 
-def test_check_consistency():
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+def test_neighborhood_search_f32():
+    neighborhood_search_test(np.float32)
+
+
+def test_neighborhood_search_f64():
+    neighborhood_search_test(np.float64)
+
+
+def check_consistency_test(dtype):
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
 
     reconstruction = pysplashsurf.reconstruct_surface(
         particles,
@@ -208,8 +224,16 @@ def test_check_consistency():
     # TODO: Delete some triangles and check for failure
 
 
-def test_tris_to_quads():
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+def test_check_consistency_f32():
+    check_consistency_test(np.float32)
+
+
+def test_check_consistency_f64():
+    check_consistency_test(np.float64)
+
+
+def tris_to_quads_test(dtype):
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
 
     mesh_with_data, reconstruction = pysplashsurf.reconstruction_pipeline(
         particles,
@@ -251,8 +275,16 @@ def test_tris_to_quads():
     assert "wnn" in mesh_with_data.point_attributes
 
 
-def test_interpolator():
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+def test_tris_to_quads_f32():
+    tris_to_quads_test(np.float32)
+
+
+def test_tris_to_quads_f64():
+    tris_to_quads_test(np.float64)
+
+
+def interpolator_test(dtype):
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
 
     mesh_with_data, reconstruction = pysplashsurf.reconstruction_pipeline(
         particles,
@@ -280,20 +312,20 @@ def test_interpolator():
     )
 
     assert type(mesh_densities) is np.ndarray
-    assert mesh_densities.dtype == np.float32
+    assert mesh_densities.dtype == dtype
     assert mesh_densities.shape == (len(mesh.vertices),)
     assert mesh_densities.min() >= 0.0
 
     mesh_particles = interpolator.interpolate_quantity(particles, mesh.vertices)
 
     assert type(mesh_particles) is np.ndarray
-    assert mesh_particles.dtype == np.float32
+    assert mesh_particles.dtype == dtype
     assert mesh_particles.shape == (len(mesh.vertices), 3)
 
     mesh_sph_normals = interpolator.interpolate_normals(mesh.vertices)
 
     assert type(mesh_sph_normals) is np.ndarray
-    assert mesh_sph_normals.dtype == np.float32
+    assert mesh_sph_normals.dtype == dtype
     assert mesh_sph_normals.shape == (len(mesh.vertices), 3)
 
     mesh_with_data.add_point_attribute("density", mesh_densities)
@@ -307,3 +339,11 @@ def test_interpolator():
     assert np.array_equal(mesh_with_data.point_attributes["density"], mesh_densities)
     assert np.array_equal(mesh_with_data.point_attributes["position"], mesh_particles)
     assert np.array_equal(mesh_with_data.point_attributes["normal"], mesh_sph_normals)
+
+
+def test_interpolator_f32():
+    interpolator_test(np.float32)
+
+
+def test_interpolator_f64():
+    interpolator_test(np.float64)

pysplashsurf/tests/test_bgeo.py

@@ -7,7 +7,15 @@
 BGEO_PATH = DIR.joinpath("ParticleData_Fluid_50.bgeo")
 
 
-def test_bgeo():
-    particles = np.array(meshio.read(BGEO_PATH).points, dtype=np.float32)
+def bgeo_test(dtype):
+    particles = np.array(meshio.read(BGEO_PATH).points, dtype=dtype)
 
     assert len(particles) == 4732
+
+
+def test_bgeo_f32():
+    bgeo_test(np.float32)
+
+
+def test_bgeo_f64():
+    bgeo_test(np.float64)

pysplashsurf/tests/test_calling.py

@@ -15,10 +15,10 @@ def now_s():
     return time.process_time_ns() / (10**9)
 
 
-def test_marching_cubes_calls():
+def marching_cubes_calls(dtype):
     print("\nTesting marching cubes calls")
 
-    particles = np.array(meshio.read(VTK_PATH).points, dtype=np.float32)
+    particles = np.array(meshio.read(VTK_PATH).points, dtype=dtype)
     reconstruction = pysplashsurf.reconstruct_surface(
         particles,
         particle_radius=0.025,
@@ -40,6 +40,14 @@ def test_marching_cubes_calls():
     assert verts_after < verts_before
 
 
+def test_marching_cubes_calls_f32():
+    marching_cubes_calls(np.float32)
+
+
+def test_marching_cubes_calls_f64():
+    marching_cubes_calls(np.float64)
+
+
 def reconstruction_pipeline(
     input_file,
     output_file,
@@ -76,9 +84,10 @@ def reconstruction_pipeline(
     quad_max_interior_angle=135.0,
     subdomain_grid=False,
     subdomain_num_cubes_per_dim=64,
+    dtype
 ):
     mesh = meshio.read(input_file)
-    particles = np.array(mesh.points, dtype=np.float64)
+    particles = np.array(mesh.points, dtype=dtype)
 
     if attributes_to_interpolate is None:
         attributes_to_interpolate = []
@@ -88,7 +97,7 @@ def reconstruction_pipeline(
     for attr in attributes_to_interpolate:
         if attr in mesh.point_data:
             if mesh.point_data[attr].dtype.kind == "f":
-                attrs[attr] = mesh.point_data[attr].astype(np.float64)
+                attrs[attr] = mesh.point_data[attr].astype(dtype)
             else:
                 attrs[attr] = mesh.point_data[attr].astype(np.int64)
 
@@ -131,11 +140,11 @@ def reconstruction_pipeline(
     mesh_with_data.write_to_file(output_file)
 
 
-def test_no_post_processing():
+def no_post_processing_test(dtype):
     start = now_s()
     subprocess.run(
         [BINARY_PATH]
-        + f"reconstruct {VTK_PATH} -o {DIR.joinpath('test_bin.vtk')} -r=0.025 -l=2.0 -c=0.5 -t=0.6 -d=on --subdomain-grid=on --mesh-cleanup=off --mesh-smoothing-weights=off --mesh-smoothing-iters=0 --normals=off --normals-smoothing-iters=0".split(),
+        + f"reconstruct {VTK_PATH} -o {DIR.joinpath('test_bin.vtk')} -r=0.025 -l=2.0 -c=0.5 -t=0.6 {"-d=on" if dtype == np.float64 else ""} --subdomain-grid=on --mesh-cleanup=off --mesh-smoothing-weights=off --mesh-smoothing-iters=0 --normals=off --normals-smoothing-iters=0".split(),
         check=True,
     )
     print("Binary done in", now_s() - start)
@@ -144,24 +153,25 @@ def test_no_post_processing():
     reconstruction_pipeline(
         VTK_PATH,
         DIR.joinpath("test.vtk"),
-        particle_radius=np.float64(0.025),
-        smoothing_length=np.float64(2.0),
-        cube_size=np.float64(0.5),
-        iso_surface_threshold=np.float64(0.6),
+        particle_radius=0.025,
+        smoothing_length=2.0,
+        cube_size=0.5,
+        iso_surface_threshold=0.6,
         mesh_smoothing_weights=False,
         mesh_smoothing_iters=0,
         normals_smoothing_iters=0,
         mesh_cleanup=False,
         compute_normals=False,
         subdomain_grid=True,
+        dtype=dtype
     )
     print("Python done in", now_s() - start)
 
     binary_mesh = meshio.read(DIR.joinpath("test_bin.vtk"))
     python_mesh = meshio.read(DIR.joinpath("test.vtk"))
 
-    binary_verts = np.array(binary_mesh.points, dtype=np.float64)
-    python_verts = np.array(python_mesh.points, dtype=np.float64)
+    binary_verts = np.array(binary_mesh.points, dtype=dtype)
+    python_verts = np.array(python_mesh.points, dtype=dtype)
 
     print("# of vertices binary:", len(binary_verts))
     print("# of vertices python:", len(python_verts))
@@ -174,11 +184,19 @@ def test_no_post_processing():
     assert np.allclose(binary_verts, python_verts)
 
 
-def test_with_post_processing():
+def test_no_post_processing_f32():
+    no_post_processing_test(np.float32)
+
+
+def test_no_post_processing_f64():
+    no_post_processing_test(np.float64)
+
+
+def with_post_processing_test(dtype):
     start = now_s()
     subprocess.run(
         [BINARY_PATH]
-        + f"reconstruct {VTK_PATH} -o {DIR.joinpath('test_bin.vtk')} -r=0.025 -l=2.0 -c=0.5 -t=0.6 -d=on --subdomain-grid=on --interpolate_attribute velocity --decimate-barnacles=on --mesh-cleanup=on --mesh-smoothing-weights=on --mesh-smoothing-iters=25 --normals=on --normals-smoothing-iters=10 --output-smoothing-weights=on --generate-quads=off".split(),
+        + f"reconstruct {VTK_PATH} -o {DIR.joinpath('test_bin.vtk')} -r=0.025 -l=2.0 -c=0.5 -t=0.6 {"-d=on" if dtype == np.float64 else ""} --subdomain-grid=on --interpolate_attribute velocity --decimate-barnacles=on --mesh-cleanup=on --mesh-smoothing-weights=on --mesh-smoothing-iters=25 --normals=on --normals-smoothing-iters=10 --output-smoothing-weights=on --generate-quads=off".split(),
         check=True,
     )
     print("Binary done in", now_s() - start)
@@ -188,12 +206,12 @@ def test_with_post_processing():
         VTK_PATH,
         DIR.joinpath("test.vtk"),
         attributes_to_interpolate=["velocity"],
-        particle_radius=np.float64(0.025),
-        smoothing_length=np.float64(2.0),
-        cube_size=np.float64(0.5),
-        iso_surface_threshold=np.float64(0.6),
+        particle_radius=0.025,
+        smoothing_length=2.0,
+        cube_size=0.5,
+        iso_surface_threshold=0.6,
         mesh_smoothing_weights=True,
-        mesh_smoothing_weights_normalization=np.float64(13.0),
+        mesh_smoothing_weights_normalization=13.0,
         mesh_smoothing_iters=25,
         normals_smoothing_iters=10,
         generate_quads=False,
@@ -203,15 +221,16 @@ def test_with_post_processing():
         decimate_barnacles=True,
         output_mesh_smoothing_weights=True,
         output_raw_normals=True,
+        dtype=dtype
     )
     print("Python done in", now_s() - start)
 
     binary_mesh = meshio.read(DIR.joinpath("test_bin.vtk"))
     python_mesh = meshio.read(DIR.joinpath("test.vtk"))
 
     # Compare number of vertices
-    binary_verts = np.array(binary_mesh.points, dtype=np.float64)
-    python_verts = np.array(python_mesh.points, dtype=np.float64)
+    binary_verts = np.array(binary_mesh.points, dtype=dtype)
+    python_verts = np.array(python_mesh.points, dtype=dtype)
 
     print("# of vertices binary:", len(binary_verts))
     print("# of vertices python:", len(python_verts))
@@ -249,3 +268,11 @@ def test_with_post_processing():
     print("Python verts:", python_verts)
 
     assert np.allclose(binary_verts, python_verts)
+
+
+def test_with_post_processing_f32():
+    with_post_processing_test(np.float32)
+
+
+def test_with_post_processing_f64():
+    with_post_processing_test(np.float64)

pysplashsurf/tests/test_sdf.py

@@ -2,7 +2,7 @@
 import numpy as np
 
 
-def test_sphere_sdf_mc():
+def sphere_sdf_mc_test(dtype):
     radius = 1.0
     num_verts = 100
 
@@ -12,7 +12,7 @@ def test_sphere_sdf_mc():
     translation = -0.5 * grid_size
 
     def make_sdf():
-        coords = np.arange(num_verts, dtype=np.float32) * dx + translation
+        coords = np.arange(num_verts, dtype=dtype) * dx + translation
         x, y, z = np.meshgrid(coords, coords, coords, indexing="ij")
         sdf = np.sqrt(x**2 + y**2 + z**2) - radius
         return sdf
@@ -31,3 +31,11 @@ def make_sdf():
     assert norms.max() < radius + 1e-4
 
     assert pysplashsurf.check_mesh_consistency(mesh, grid) is None
+
+
+def test_sphere_sdf_mc_f32():
+    sphere_sdf_mc_test(np.float32)
+
+
+def test_sphere_sdf_mc_f64():
+    sphere_sdf_mc_test(np.float64)