|
1 | | -#  |
2 | | -[](https://crates.io/crates/splashsurf) |
3 | | -[](https://docs.rs/splashsurf_lib) |
4 | | -[](https://github.com/w1th0utnam3/splashsurf/blob/master/LICENSE) |
5 | | - |
| 1 | +# splashsurf_lib |
| 2 | +Library for surface reconstruction of SPH particle data |
6 | 3 |
|
7 | | -Surface reconstruction library and CLI for particle data from SPH simulations, written in Rust. |
| 4 | +The library is mainly used by the `splashsurf` command-line tool which is [also available](https://crates.io/crates/splashsurf) on crates.io. |
| 5 | + |
| 6 | +For more information about the CLI, check out the [readme in the root of the repository](https://github.com/w1th0utnam3/splashsurf). |
8 | 7 |
|
9 | 8 | **Contents** |
10 | | -- [!splashsurf logo](#) |
11 | | -- [The `splashsurf` CLI](#the-splashsurf-cli) |
12 | | - - [Introduction](#introduction) |
13 | | - - [Notes](#notes) |
14 | | - - [Installation](#installation) |
| 9 | +- [splashsurf_lib](#splashsurf_lib) |
15 | 10 | - [Usage](#usage) |
16 | | - - [Basic usage](#basic-usage) |
17 | | - - [Sequences of files](#sequences-of-files) |
18 | | - - [Input file formats](#input-file-formats) |
19 | | - - [VTK](#vtk) |
20 | | - - [PLY](#ply) |
21 | | - - [XYZ](#xyz) |
22 | | - - [Output file formats](#output-file-formats) |
23 | | -- [License](#license) |
24 | | - |
25 | | -# The `splashsurf` CLI |
26 | | - |
27 | | -The following sections mainly focus on the CLI of `splashsurf`. For more information on the library parts, see the [corresponding readme](https://github.com/w1th0utnam3/splashsurf/blob/master/splashsurf_lib/README.md) in the `splashsurf_lib` subfolder. |
28 | | - |
29 | | -## Introduction |
30 | | - |
31 | | -This is a basic but high-performance implementation of a marching cubes based surface reconstruction for SPH fluid simulations (e.g performed with [SPlisHSPlasH](https://github.com/InteractiveComputerGraphics/SPlisHSPlasH)). |
32 | | -The output of this tool is the reconstructed triangle surface mesh of the fluid. |
33 | | -At the moment it does not compute normals or other additional data. |
34 | | -As input, it supports particle positions from .vtk files and binary .xyz files (i.e. files containing a binary dump of a particle position array). In addition, required parameters are the kernel radius and particle radius (to compute the volume of particles) used for the original SPH simulation as well as the surface threshold. |
35 | | - |
36 | | -The implementation first computes the density of each particle using the typical SPH approach with a cubic kernel. |
37 | | -This density is then evaluated or mapped onto a sparse grid using spatial hashing in the support radius of each particle. |
38 | | -This implies that memory is only allocated in areas where the fluid density is non-zero. This is in contrast to a naive approach where the marching cubes background grid is allocated for the whole domain. |
39 | | -Finally, the marching cubes reconstruction is performed only in those grid cells where the density values cross the surface threshold. Cells completely in the interior of the fluid are skipped. For more details, please refer to the [readme of the library]((https://github.com/w1th0utnam3/splashsurf/blob/master/splashsurf_lib/README.md)). |
40 | | - |
41 | | -## Notes |
42 | | - |
43 | | -Due the use of hash maps and multi-threading (if enabled), the output of this implementation is not deterministic. |
44 | | -In the future, flags may be added to switch the internal data structures to use binary trees for debugging purposes. |
45 | | - |
46 | | -Note that for small numbers of fluid particles (i.e. in the low thousands or less) the multi-threaded implementation may have worse performance due to the worker pool overhead and looks on the map data structures (even though the library uses [`dashmap`](https://crates.io/crates/dashmap) which is more optimized for multi-threaded usage). |
47 | | - |
48 | | -As shown below, the tool can handle the output of large simulations. |
49 | | -However, it was not tested with a wide range of parameters and may not be totally robust against corner-cases or extreme parameters. |
50 | | -If you experience problems, please report them together with your input data. |
51 | | - |
52 | | -## Installation |
53 | | - |
54 | | -The command-line tool can be built from this repository but is also available on crates.io. |
55 | | -If you have a [Rust toolchain installed](https://www.rust-lang.org/tools/install) you can install `splashsurf` with the command |
56 | | -``` |
57 | | -cargo install splashsurf |
58 | | -``` |
| 11 | + - [Feature flags](#feature-flags) |
| 12 | + - [The surface reconstruction procedure](#the-surface-reconstruction-procedure) |
59 | 13 |
|
60 | 14 | ## Usage |
61 | 15 |
|
62 | | -### Basic usage |
63 | | - |
64 | | -``` |
65 | | -USAGE: |
66 | | - splashsurf.exe [FLAGS] [OPTIONS] <input-file> --cube-size <cube-size> --kernel-radius <kernel-radius> --particle-radius <particle-radius> --surface-threshold <surface-threshold> |
67 | | -
|
68 | | -FLAGS: |
69 | | - -h, --help Prints help information |
70 | | - --mt-files Whether to enable multi-threading to process multiple input files in parallel, conflicts with |
71 | | - --mt-particles |
72 | | - --mt-particles Whether to enable multi-threading for a single input file by processing chunks of particles in |
73 | | - parallel, conflicts with --mt-files |
74 | | - -d Whether to enable the use of double precision for all computations (disabled by default) |
75 | | - -V, --version Prints version information |
76 | | -
|
77 | | -OPTIONS: |
78 | | - --cube-size <cube-size> |
79 | | - The marching cubes grid size in multiplies of the particle radius |
80 | | -
|
81 | | - --domain-max <domain-max> <domain-max> <domain-max> |
82 | | - Upper corner of the domain where surface reconstruction should be performed, format:domain- |
83 | | - max=x_max;y_max;z_max (requires domain-min to be specified) |
84 | | - --domain-min <domain-min> <domain-min> <domain-min> |
85 | | - Lower corner of the domain where surface reconstruction should be performed, format: domain- |
86 | | - min=x_min;y_min;z_min (requires domain-max to be specified) |
87 | | - --kernel-radius <kernel-radius> |
88 | | - The kernel radius used for building the density map in multiplies of the particle radius |
89 | | -
|
90 | | - --output-dir <output-dir> Optional base directory for all output files |
91 | | - --output-dm-grid <output-dm-grid> |
92 | | - Optional filename for writing the grid representation of the intermediate density map to disk |
93 | | -
|
94 | | - --output-dm-points <output-dm-points> |
95 | | - Optional filename for writing the point cloud representation of the intermediate density map to disk |
96 | | -
|
97 | | - -o <output-file> Filename for writing the reconstructed surface to disk |
98 | | - --particle-radius <particle-radius> The particle radius of the input data |
99 | | - --rest-density <rest-density> The rest density of the fluid [default: 1000.0] |
100 | | - --splash-detection-radius <splash-detection-radius> |
101 | | - If a particle has no neighbors in this radius (in multiplies of the particle radius) it is considered as a |
102 | | - free particles |
103 | | - --surface-threshold <surface-threshold> |
104 | | - The iso-surface threshold for the density, i.e. value of the reconstructed density that indicates the fluid |
105 | | - surface |
106 | | -
|
107 | | -ARGS: |
108 | | - <input-file> Path to the input file where the particle positions are stored (supported formats: VTK, binary |
109 | | - XYZ, PLY) |
110 | | -``` |
111 | | -For example: |
112 | | -``` |
113 | | -splashsurf "/home/user/canyon.xyz" --output-dir="/home/user/temp" --particle-radius=0.011 --kernel-radius=4.0 --cube-size=1.5 --surface-threshold=0.6 --mt-particles |
114 | | -``` |
115 | | -With these parameters, a scene with 13353401 particles is reconstructed in nearly than 25s on a i9 9900K. The output is a mesh with 6016212 triangles. |
116 | | -``` |
117 | | -[2020-08-25T15:52:34.515885+02:00][splashsurf::reconstruction][INFO] Loading dataset from "/local/data/temp/canyon.xyz"... |
118 | | -[2020-08-25T15:52:34.655354+02:00][splashsurf::reconstruction][INFO] Loaded dataset with 13353401 particle positions. |
119 | | -[2020-08-25T15:52:34.684734+02:00][splashsurf_lib][INFO] Minimal enclosing bounding box of particles was computed as: AxisAlignedBoundingBox { min: [-25.0060978, -5.0146289, -40.0634613], max: [24.4994926, 18.3062096, 39.7757950] } |
120 | | -[2020-08-25T15:52:34.684748+02:00][splashsurf_lib][INFO] Using a grid with 6002x2828x9679 points and 6001x2827x9678 cells of edge length 0.0165. |
121 | | -[2020-08-25T15:52:34.684753+02:00][splashsurf_lib][INFO] The resulting domain size is: AxisAlignedBoundingBox { min: [-49.7588959, -16.6750488, -79.9830933], max: [49.2576065, 29.9704514, 79.7039032] } |
122 | | -[2020-08-25T15:52:34.684756+02:00][splashsurf_lib][INFO] Starting neighborhood search... |
123 | | -[2020-08-25T15:52:36.570860+02:00][splashsurf_lib][INFO] Computing particle densities... |
124 | | -[2020-08-25T15:52:37.645919+02:00][splashsurf_lib::density_map][INFO] Starting construction of sparse density map for 13353401 particles... |
125 | | -[2020-08-25T15:52:37.653068+02:00][splashsurf_lib::density_map][INFO] To take into account the kernel evaluation radius, the allowed domain of particles was restricted to: AxisAlignedBoundingBox { min: [-49.7093468, -16.6254997, -79.9335403], max: [49.2080574, 29.9209023, 79.6543503] } |
126 | | -[2020-08-25T15:52:55.559939+02:00][splashsurf_lib::density_map][INFO] Sparse density map with 31519986 point data values was constructed. |
127 | | -[2020-08-25T15:52:55.559999+02:00][splashsurf_lib::density_map][INFO] Construction of sparse density map done. |
128 | | -[2020-08-25T15:52:55.560005+02:00][splashsurf_lib::marching_cubes][INFO] Starting interpolation of cell data for marching cubes... |
129 | | -[2020-08-25T15:52:59.118442+02:00][splashsurf_lib::marching_cubes][INFO] Generated cell data for marching cubes with 3009863 cells and 3011516 vertices. |
130 | | -[2020-08-25T15:52:59.118470+02:00][splashsurf_lib::marching_cubes][INFO] Interpolation done. |
131 | | -[2020-08-25T15:52:59.118474+02:00][splashsurf_lib::marching_cubes][INFO] Starting marching cubes triangulation of 3009863 cells... |
132 | | -[2020-08-25T15:52:59.279570+02:00][splashsurf_lib::marching_cubes][INFO] Generated surface mesh with 6016212 triangles and 3011516 vertices. |
133 | | -[2020-08-25T15:52:59.279597+02:00][splashsurf_lib::marching_cubes][INFO] Triangulation done. |
134 | | -[2020-08-25T15:52:59.296979+02:00][splashsurf::reconstruction][INFO] Writing surface mesh to "/home/floeschner/programming/temp/canyon_surface.xyz"... |
135 | | -[2020-08-25T15:52:59.808101+02:00][splashsurf::reconstruction][INFO] Done. |
136 | | -[2020-08-25T15:52:59.879069+02:00][splashsurf][INFO] Finished processing all inputs. |
137 | | -[2020-08-25T15:52:59.879103+02:00][splashsurf][INFO] surface reconstruction cli: 100.00%, 25363.19ms/call @ 0.04Hz |
138 | | -[2020-08-25T15:52:59.879107+02:00][splashsurf][INFO] loading particle positions: 0.55%, 139.30ms/call @ 0.04Hz |
139 | | -[2020-08-25T15:52:59.879109+02:00][splashsurf][INFO] reconstruct_surface: 97.15%, 24641.60ms/call @ 0.04Hz |
140 | | -[2020-08-25T15:52:59.879111+02:00][splashsurf][INFO] compute minimum enclosing aabb: 0.12%, 29.37ms/call @ 0.04Hz |
141 | | -[2020-08-25T15:52:59.879113+02:00][splashsurf][INFO] neighborhood_search: 7.65%, 1886.10ms/call @ 0.04Hz |
142 | | -[2020-08-25T15:52:59.879115+02:00][splashsurf][INFO] parallel_generate_cell_to_particle_map: 11.52%, 217.26ms/call @ 0.04Hz |
143 | | -[2020-08-25T15:52:59.879117+02:00][splashsurf][INFO] get_cell_neighborhoods_par: 2.70%, 50.97ms/call @ 0.04Hz |
144 | | -[2020-08-25T15:52:59.879119+02:00][splashsurf][INFO] calculate_particle_neighbors_par: 52.80%, 995.77ms/call @ 0.04Hz |
145 | | -[2020-08-25T15:52:59.879121+02:00][splashsurf][INFO] parallel_compute_particle_densities: 1.12%, 275.39ms/call @ 0.04Hz |
146 | | -[2020-08-25T15:52:59.879144+02:00][splashsurf][INFO] parallel_generate_sparse_density_map: 72.70%, 17914.09ms/call @ 0.04Hz |
147 | | -[2020-08-25T15:52:59.879146+02:00][splashsurf][INFO] triangulate_density_map: 15.17%, 3736.97ms/call @ 0.04Hz |
148 | | -[2020-08-25T15:52:59.879148+02:00][splashsurf][INFO] interpolate_points_to_cell_data: 95.22%, 3558.47ms/call @ 0.04Hz |
149 | | -[2020-08-25T15:52:59.879161+02:00][splashsurf][INFO] generate_iso_surface_vertices: 69.82%, 2484.44ms/call @ 0.04Hz |
150 | | -[2020-08-25T15:52:59.879164+02:00][splashsurf][INFO] relative_to_threshold_postprocessing: 29.70%, 1056.71ms/call @ 0.04Hz |
151 | | -[2020-08-25T15:52:59.879167+02:00][splashsurf][INFO] triangulate: 4.78%, 178.50ms/call @ 0.04Hz |
152 | | -[2020-08-25T15:52:59.879171+02:00][splashsurf][INFO] write surface mesh to file: 2.02%, 511.21ms/call @ 0.04Hz |
153 | | -``` |
154 | | - |
155 | | -### Sequences of files |
156 | | - |
157 | | -*TODO: Describe syntax to reconstruct a sequence of files* |
158 | | - |
159 | | -## Input file formats |
160 | | - |
161 | | -### VTK |
162 | | - |
163 | | -Files with the "`.vtk`" extension are loaded using [`vtkio`](https://crates.io/crates/vtkio). The VTK file is loaded as a big endian binary file and has to contain an "Unstructured Grid" with either `f32` or `f64` vertex coordinates. Any other data or attributes are ignored. Only the first "Unstructured Grid" is loaded, other entities are ignored. |
164 | | - |
165 | | -### PLY |
166 | | - |
167 | | -Files with the "`.ply`" extension are loaded using [`ply-rs`](https://crates.io/crates/ply-rs). The PLY file has to contain an element called "`vertex`" with the properties `x`, `y` and `z` of type `f32`/["`Property::Float`"](https://docs.rs/ply-rs/0.1.3/ply_rs/ply/enum.Property.html#variant.Float). Any other properties or elements are ignored. |
168 | | - |
169 | | -### XYZ |
170 | | - |
171 | | -Files with the "`.xyz`" extension are interpreted as raw bytes of `f32` values in native endianness of the system. Three consecutive `f32`s represent a (x,y,z) coordinate triplet of a fluid particle. |
172 | | - |
173 | | -## Output file formats |
174 | | - |
175 | | -Currently, only VTK files are supported for output. |
176 | | - |
177 | | -# License |
178 | | - |
179 | | -For license information of this project, see the LICENSE file. |
180 | | -The splashsurf logo is based on two graphics ([1](https://www.svgrepo.com/svg/295647/wave), [2](https://www.svgrepo.com/svg/295652/surfboard-surfboard)) published on SVG Repo under a CC0 ("No Rights Reserved") license. |
| 16 | +A reconstruction from particle positions is performed using the top-level `reconstruct_surface` function: |
| 17 | +```rust |
| 18 | +pub fn reconstruct_surface<I: Index, R: Real>( |
| 19 | + particle_positions: &[Vector3<R>], |
| 20 | + parameters: &Parameters<R>, |
| 21 | +) -> Result<SurfaceReconstruction<I, R>, ReconstructionError<I, R>> |
| 22 | +``` |
| 23 | +See the [documentation of the crate](https://docs.rs/splashsurf_lib/latest/splashsurf_lib/) on docs.rs for more information on the usage. |
| 24 | + |
| 25 | +The library re-exports `nalgebra` to avoid version conflicts for users of the library. |
| 26 | + |
| 27 | +## Feature flags |
| 28 | + |
| 29 | +By default none of the following features are enabled to reduce the dependencies introduced by this library. The following feature flags are available for `splashsurf_lib`: |
| 30 | + |
| 31 | + - **vtk-extras**: Enables convenience traits and helper functions to convert the mesh types returned by the library to [`vtkio`](https://crates.io/crates/vtkio) data structures (in particular [`UnstructuredGridPiece`](https://docs.rs/vtkio/latest/vtkio/model/struct.UnstructuredGridPiece.html)) that can be used to easily write the meshes to VTK files (e.g. for viewing them with [Paraview](https://www.paraview.org/)). Check out the documentation of `vtkio` or the [corresponding io module](https://github.com/w1th0utnam3/splashsurf/blob/master/splashsurf/src/io/vtk_format.rs) of the `splashsurf` CLI for reference. |
| 32 | + - **profiling**: Enables profiling of the library using [`coarse-prof`](https://crates.io/crates/coarse-prof). Several functions in the library will use the [`profile!`](https://docs.rs/coarse-prof/latest/coarse_prof/macro.profile.html) macro with the function name as an argument to record their runtime. The user of the library can then obtain the profiling data using the functions provided by the `coarse-prof` crate. Note that profiling using this crate might reduce performance for surface reconstructions with a very small number of particles (i.e. only a few hundred). |
| 33 | + |
| 34 | +For each of the features, `splashsurf_lib` re-exports the corresponding dependencies to avoid version conflicts for users of the library. |
| 35 | + |
| 36 | +## The surface reconstruction procedure |
| 37 | + |
| 38 | +Currently, only one method based on a "spatial hashing" strategy is implemented. |
| 39 | + |
| 40 | +**Short summary**: The fluid density is evaluated or mapped onto a sparse grid using spatial hashing in the support radius of each fluid particle. This implies that memory is only allocated in areas where the fluid density is non-zero. This is in contrast to a naive approach where the marching cubes background grid is allocated for the whole domain. Finally, the marching cubes reconstruction is performed only in those grid cells where an edge crosses the surface threshold. Cells completely in the interior of the fluid are skipped in the marching cubes phase. |
| 41 | + |
| 42 | +**Individual steps**: |
| 43 | + 1. Construct a "virtual background grid" with the desired resolution of the marching cubes algorithm. In the end, the procedure will place a single surface mesh vertex on each edge of this virtual grid, where the fluid surface crosses the edge (or rather, where the fluid density crosses the specified threshold). Virtual means that no storage is actually allocated for this grid yet; only its topology is used implicitly later. |
| 44 | + 2. Compute the density of each fluid particle |
| 45 | + - Perform a neighborhood search |
| 46 | + - Per particle, evaluate an SPH sum over the neighbors to compute its density (based on input parameters of kernel radius and particle rest mass) |
| 47 | + 3. Optional: filter out (or rather mask as inactive) single particles if the user provided a "splash detection radius". This is done by performing an additional neighborhood search using this splash detection radius instead of the kernel radius. |
| 48 | + 4. Compute a "sparse density map": a map from the index of a vertex of the virtual background grid to the corresponding fluid density value. The map will only contain entries for vertices where the fluid density is non-zero. Construction of the map: |
| 49 | + - Iterate over all active particles |
| 50 | + - For each particle evaluate its kernel at all virtual background vertices that can be influenced by it (i.e. vertices inside its kernel radius) |
| 51 | + - Add-assign the corresponding density contribution (kernel value times particle density) to the vertex entry in the density map |
| 52 | + 5. Interpolate the density values at the vertices of each virtual background cell to points on the edges where the edge crosses the fluid surface |
| 53 | + - Iterate over all vertices in the sparse density map |
| 54 | + - Skip vertices where the density value is above the threshold to be considered inside of the fluid |
| 55 | + - For each of the remaining vertices, check if any neighboring vertex (i.e. a vertex that is directly connected with it by an edge of the virtual background grid) is above the threshold |
| 56 | + - If this is the case, an edge that crosses the fluid surface was found and the position of the surface on the edge can be calculated using linear interpolation |
| 57 | + - This interpolated position is stored in a new map that maps indices of virtual background cells to structs containing the surface crossing points per edge. Entries are only created for cells with edges that actually cross the surface |
| 58 | + 6. Triangulate the points on the edges using a marching cubes case table |
| 59 | + - Iterate over all cells in the cell data map |
| 60 | + - For each cell, look up the corresponding triangulation in a marching cubes LUT |
| 61 | + - Emit the required triangles into the final mesh data structure |
0 commit comments