diff --git a/.gitignore b/.gitignore index 35e928c..25a5911 100644 --- a/.gitignore +++ b/.gitignore @@ -15,3 +15,4 @@ STYLE.md Cargo.lock .idea .vscode +.DS_Store \ No newline at end of file diff --git a/Cargo.toml b/Cargo.toml index 9290235..0500a19 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -1,5 +1,6 @@ [workspace] -members = ["build/salva2d", "build/salva3d", "examples2d", "examples3d"] +members = ["build/salva2d", "build/salva2d-f64", "build/salva3d", "build/salva3d-f64", "examples2d", "examples3d"] +default-members = ["build/salva2d", "build/salva2d-f64", "build/salva3d", "build/salva3d-f64"] resolver = "2" [workspace.lints] @@ -17,7 +18,9 @@ opt-level = 1 [patch.crates-io] salva2d = { path = "./build/salva2d" } +salva2d-f64 = { path = "./build/salva2d-f64" } salva3d = { path = "./build/salva3d" } +salva3d-f64 = { path = "./build/salva3d-f64" } #parry2d = { git = "https://github.com/sebcrozet/parry" } #parry3d = { git = "https://github.com/sebcrozet/parry" } diff --git a/README.md b/README.md index 99343ab..1ab2421 100644 --- a/README.md +++ b/README.md @@ -46,3 +46,10 @@ is inspired from its renown painting [The Persistence of Memory](https://en.wiki - **Multiphase fluids**: mix several fluids with different characteristics (densities, viscosities, etc.) - Optional **two-way coupling** with bodies from **rapier**. - **WASM** support + +## Running examples + +```bash +cargo run --release -p examples2d --bin all_examples2 +cargo run --release -p examples3d --bin all_examples3 +``` diff --git a/build/salva2d-f64/Cargo.toml b/build/salva2d-f64/Cargo.toml new file mode 100644 index 0000000..8f43212 --- /dev/null +++ b/build/salva2d-f64/Cargo.toml @@ -0,0 +1,64 @@ +[package] +name = "salva2d-f64" +version = "0.9.0" +authors = [ "Sébastien Crozet " ] +description = "2-dimensional particle-based fluid dynamics in Rust." +documentation = "https://salva.rs/docs" +homepage = "https://salva.rs" +repository = "https://github.com/dimforge/salva" +readme = "README.md" +categories = [ "science", "game-development", "mathematics", "simulation", "wasm"] +keywords = [ "physics", "dynamics", "particles", "fluids", "SPH" ] +license = "Apache-2.0" +edition = "2021" + +[badges] +maintenance = { status = "actively-developed" } + +[features] +default = [ "dim2", "f64" ] +dim2 = [ ] +f64 = [ ] +parallel = [ "rayon" ] +sampling = [ "rapier" ] +rapier = [ "parry", "rapier2d-f64" ] +rapier-testbed = [ "rapier", "rapier_testbed2d", "graphics" ] +rapier-harness = [ "rapier-testbed" ] +parry = [ "parry2d-f64" ] +wasm-bindgen = [ ] +graphics = [ "bevy", "bevy_egui" ] + +[lints] +rust.unexpected_cfgs = { level = "warn", check-cfg = [ + 'cfg(feature, values("dim3", "f32"))', +] } + +[lib] +name = "salva2d_f64" +path = "../../src/lib.rs" +required-features = [ "dim2", "f64" ] + +[dependencies] +approx = "0.5" +num-traits = "0.2" +fnv = "1.0" +itertools = "0.13" +generational-arena = "0.2" +rayon = { version = "1.8", optional = true } + +nalgebra = { version = "0.35", features = ["convert-glam033"] } +glamx = "0.3" +parry2d-f64 = { version = "0.28", optional = true } +rapier2d-f64 = { version = "0.33", optional = true } +# TODO update it to f64 +rapier_testbed2d = { version = "0.33", optional = true } + +bevy_egui = { version = "0.26", features = ["immutable_ctx"], optional = true } +bitflags = "2" + +[target.'cfg(not(target_arch = "wasm32"))'.dependencies] +bevy = { version = "0.13.2", default-features = false, features = ["bevy_winit", "bevy_render", "x11"], optional = true } + +# Dependencies for WASM only. +[target.'cfg(target_arch = "wasm32")'.dependencies] +bevy = { version = "0.13", default-features = false, features = ["bevy_winit", "bevy_render"], optional = true } diff --git a/build/salva2d/Cargo.toml b/build/salva2d/Cargo.toml index b44caaf..61db660 100644 --- a/build/salva2d/Cargo.toml +++ b/build/salva2d/Cargo.toml @@ -22,23 +22,26 @@ edition = "2021" maintenance = { status = "actively-developed" } [features] -default = ["dim2"] -dim2 = [] -parallel = ["rayon"] -sampling = ["rapier"] -rapier = ["parry", "rapier2d"] -rapier-testbed = ["rapier", "rapier_testbed2d", "graphics"] -rapier-harness = ["rapier-testbed"] -parry = ["parry2d"] -graphics = ["bevy", "bevy_egui"] +default = [ "dim2", "f32" ] +dim2 = [ ] +f32 = [ ] +parallel = [ "rayon" ] +sampling = [ "rapier" ] +rapier = [ "parry", "rapier2d" ] +rapier-testbed = [ "rapier", "rapier_testbed2d", "graphics", "kiss3d" ] +rapier-harness = [ "rapier-testbed" ] +parry = [ "parry2d" ] +graphics = [ "bevy", "bevy_egui" ] [lints] -workspace = true +rust.unexpected_cfgs = { level = "warn", check-cfg = [ + 'cfg(feature, values("dim3", "f64"))', +] } [lib] name = "salva2d" path = "../../src/lib.rs" -required-features = ["dim2"] +required-features = [ "dim2", "f32" ] [dependencies] approx = "0.5" @@ -46,13 +49,14 @@ num-traits = "0.2" fnv = "1.0" itertools = "0.14" generational-arena = "0.2" -instant = { version = "0.1", features = ["now"] } rayon = { version = "1.8", optional = true } -nalgebra = "0.33" -parry2d = { version = "0.18", optional = true } -rapier2d = { version = "0.23", optional = true } -rapier_testbed2d = { version = "0.23", optional = true } +nalgebra = { version = "0.35", features = ["convert-glam033"] } +glamx = "0.3" +parry2d = { version = "0.28", optional = true } +rapier2d = { version = "0.33", optional = true } +rapier_testbed2d = { version = "0.33", optional = true } +kiss3d = { version = "0.43", optional = true } bevy_egui = { version = "0.31", features = ["immutable_ctx"], optional = true } bitflags = "2" diff --git a/build/salva3d-f64/Cargo.toml b/build/salva3d-f64/Cargo.toml new file mode 100644 index 0000000..7cfc7c6 --- /dev/null +++ b/build/salva3d-f64/Cargo.toml @@ -0,0 +1,60 @@ +[package] +name = "salva3d-f64" +version = "0.9.0" +authors = [ "Sébastien Crozet " ] +description = "3-dimensional particle-based fluid dynamics in Rust." +documentation = "https://salva.rs/rustdoc/salva3d/index.html" +homepage = "https://salva.rs" +repository = "https://github.com/dimforge/salva" +readme = "README.md" +keywords = [ "physics", "dynamics", "particles", "fluids", "SPH" ] +license = "Apache-2.0" +edition = "2021" + +[features] +default = [ "dim3", "f64" ] +dim3 = [ ] +f64 = [ ] +parallel = [ "rayon" ] +rapier = [ "parry", "rapier3d-f64" ] +sampling = [ "rapier" ] +rapier-testbed = [ "rapier", "rapier_testbed3d", "graphics" ] +rapier-harness = [ "rapier-testbed" ] +parry = [ "parry3d-f64" ] +wasm-bindgen = [ ] +graphics = [ "bevy", "bevy_egui" ] + +[lints] +rust.unexpected_cfgs = { level = "warn", check-cfg = [ + 'cfg(feature, values("dim2", "f32"))', +] } + +[lib] +name = "salva3d_f64" +path = "../../src/lib.rs" +required-features = [ "dim3", "f64" ] + +[dependencies] +approx = "0.5" +num-traits = "0.2" +fnv = "1.0" +itertools = "0.13" +generational-arena = "0.2" +rayon = { version = "1.8", optional = true } + +nalgebra = { version = "0.35", features = ["convert-glam033"] } +glamx = "0.3" +parry3d-f64 = { version = "0.28", optional = true } +rapier3d-f64 = { version = "0.33", optional = true } +# TODO update it to f64 +rapier_testbed3d = { version = "0.33", optional = true } + +bevy_egui = { version = "0.26", features = ["immutable_ctx"], optional = true } +bitflags = "2" + +[target.'cfg(not(target_arch = "wasm32"))'.dependencies] +bevy = { version = "0.13", default-features = false, features = ["bevy_winit", "bevy_render", "x11"], optional = true } + +# Dependencies for WASM only. +[target.'cfg(target_arch = "wasm32")'.dependencies] +bevy = { version = "0.13", default-features = false, features = ["bevy_winit", "bevy_render"], optional = true } diff --git a/build/salva3d/Cargo.toml b/build/salva3d/Cargo.toml index 8fe1fee..978619d 100644 --- a/build/salva3d/Cargo.toml +++ b/build/salva3d/Cargo.toml @@ -12,23 +12,26 @@ license = "Apache-2.0" edition = "2021" [lints] -workspace = true +rust.unexpected_cfgs = { level = "warn", check-cfg = [ + 'cfg(feature, values("dim2", "f64"))', +] } [features] -default = ["dim3"] -dim3 = [] -parallel = ["rayon"] -rapier = ["parry", "rapier3d"] -sampling = ["rapier"] -rapier-testbed = ["rapier", "rapier_testbed3d", "graphics"] -rapier-harness = ["rapier-testbed"] -parry = ["parry3d"] -graphics = ["bevy", "bevy_egui"] +default = [ "dim3", "f32" ] +dim3 = [ ] +f32 = [ ] +parallel = [ "rayon" ] +rapier = [ "parry", "rapier3d" ] +sampling = [ "rapier" ] +rapier-testbed = [ "rapier", "rapier_testbed3d", "graphics", "kiss3d" ] +rapier-harness = [ "rapier-testbed" ] +parry = [ "parry3d" ] +graphics = [ "bevy", "bevy_egui" ] [lib] name = "salva3d" path = "../../src/lib.rs" -required-features = ["dim3"] +required-features = [ "dim3", "f32" ] [dependencies] approx = "0.5" @@ -36,13 +39,14 @@ num-traits = "0.2" fnv = "1.0" itertools = "0.14" generational-arena = "0.2" -instant = { version = "0.1", features = ["now"] } rayon = { version = "1.8", optional = true } -nalgebra = "0.33" -parry3d = { version = "0.18", optional = true } -rapier3d = { version = "0.23", optional = true } -rapier_testbed3d = { version = "0.23.1", optional = true } +nalgebra = { version = "0.35", features = ["convert-glam033"] } +glamx = "0.3" +parry3d = { version = "0.28", optional = true } +rapier3d = { version = "0.33", optional = true } +rapier_testbed3d = { version = "0.33", optional = true } +kiss3d = { version = "0.43", optional = true } bevy_egui = { version = "0.31", features = ["immutable_ctx"], optional = true } bitflags = "2.6.0" diff --git a/examples2d/Cargo.toml b/examples2d/Cargo.toml index 4bfe76f..49c1820 100644 --- a/examples2d/Cargo.toml +++ b/examples2d/Cargo.toml @@ -5,17 +5,18 @@ authors = ["Sébastien Crozet "] edition = "2018" [features] -default = [] -parallel = ["rapier_testbed2d/parallel"] +default = ["parallel"] +parallel = [ "rapier_testbed2d/parallel", "salva2d/parallel" ] [dependencies] -Inflector = "0.11" -nalgebra = "0.33" -parry2d = "0.18" -rapier2d = "0.23" -rapier_testbed2d = "0.23" -parry3d = "0.18" +Inflector = "0.11" +nalgebra = "0" +parry2d = "0.28" +rapier2d = "0.33" +rapier_testbed2d = "0.33" +parry3d = "0.28" bevy = "0.15" +pollster = "0.4" [dependencies.salva2d] path = "../build/salva2d" diff --git a/examples2d/all_examples2.rs b/examples2d/all_examples2.rs index 4b93d43..7d69dfb 100644 --- a/examples2d/all_examples2.rs +++ b/examples2d/all_examples2.rs @@ -4,7 +4,7 @@ extern crate nalgebra as na; use inflector::Inflector; -use rapier_testbed2d::{Testbed, TestbedApp}; +use rapier_testbed2d::{Example, TestbedApp}; mod basic2; mod custom_forces2; @@ -46,20 +46,21 @@ fn main() { .unwrap_or(String::new()) .to_camel_case(); - let mut builders: Vec<(_, fn(&mut Testbed))> = vec![ - ("Basic", basic2::init_world), - ("Layers", layers2::init_world), - ("Custom forces", custom_forces2::init_world), - ("Elasticity", elasticity2::init_world), - ("Surface tension", surface_tension2::init_world), + let mut builders = vec![ + Example::demo("Basic", basic2::init_world), + Example::demo("Layers", layers2::init_world), + Example::demo("Custom forces", custom_forces2::init_world), + Example::demo("Elasticity", elasticity2::init_world), + Example::demo("Surface tension", surface_tension2::init_world), ]; - builders.sort_by_key(|builder| builder.0); + builders.sort_by_key(|builder| builder.name); let i = builders .iter() - .position(|builder| builder.0.to_camel_case().as_str() == demo.as_str()) + .position(|builder| builder.name.to_camel_case().as_str() == demo.as_str()) .unwrap_or(0); - let testbed = TestbedApp::from_builders(i, builders); + builders.rotate_left(i); + let testbed = TestbedApp::from_builders(builders); - testbed.run() + pollster::block_on(testbed.run()); } diff --git a/examples2d/basic2.rs b/examples2d/basic2.rs index cda0203..e077f5a 100644 --- a/examples2d/basic2.rs +++ b/examples2d/basic2.rs @@ -1,6 +1,6 @@ extern crate nalgebra as na; -use na::{DVector, Point2, Point3, Vector2}; +use na::{Vector2, Vector3}; use rapier2d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier2d::geometry::{Collider, ColliderBuilder, ColliderSet}; use rapier_testbed2d::Testbed; @@ -38,8 +38,8 @@ pub fn init_world(testbed: &mut Testbed) { for j in 0..nj { let x = (i as f32) * PARTICLE_RADIUS * 2.0 - ni as f32 * PARTICLE_RADIUS; let y = (j as f32 + 1.0) * PARTICLE_RADIUS * 2.0 + 0.5; - points1.push(Point2::new(x, y)); - points2.push(Point2::new(x + ni as f32 * PARTICLE_RADIUS, y)); + points1.push(Vector2::new(x, y)); + points2.push(Vector2::new(x + ni as f32 * PARTICLE_RADIUS, y)); } } @@ -47,7 +47,7 @@ pub fn init_world(testbed: &mut Testbed) { for j in 0..nj * 2 { let x = (i as f32) * PARTICLE_RADIUS * 2.0 - ni as f32 * PARTICLE_RADIUS; let y = (j as f32 + 1.0) * PARTICLE_RADIUS * 2.0 + 0.5; - points3.push(Point2::new(x, y + shift2)); + points3.push(Vector2::new(x, y + shift2)); } } @@ -57,7 +57,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); let elasticity: Becker2009Elasticity = Becker2009Elasticity::new(1_000.0, 0.3, true); let viscosity = XSPHViscosity::new(0.5, 1.0); @@ -65,13 +65,13 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(1.0, 0.4, 0.6)); + plugin.set_fluid_color(fluid_handle, Vector3::new(1.0, 0.4, 0.6)); let viscosity = ArtificialViscosity::new(0.5, 0.0); let mut fluid = Fluid::new(points3, PARTICLE_RADIUS, 1.0, InteractionGroups::default()); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.6, 0.8, 0.5)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.6, 0.8, 0.5)); /* * Ground @@ -79,17 +79,19 @@ pub fn init_world(testbed: &mut Testbed) { let ground_size = Vector2::new(10.0, 1.0); let nsubdivs = 50; - let heights = DVector::from_fn(nsubdivs + 1, |i, _| { - if i == 0 || i == nsubdivs { - 20.0 - } else { - (i as f32 * ground_size.x / (nsubdivs as f32)).cos() * 0.5 - } - }); + let heights: Vec<_> = (0..=nsubdivs) + .map(|i| { + if i == 0 || i == nsubdivs { + 20.0 + } else { + (i as f32 * ground_size.x / (nsubdivs as f32)).cos() * 0.5 + } + }) + .collect(); let rigid_body = RigidBodyBuilder::fixed().build(); let handle = bodies.insert(rigid_body); - let collider = ColliderBuilder::heightfield(heights, ground_size).build(); + let collider = ColliderBuilder::heightfield(heights, ground_size.into()).build(); let co_handle = colliders.insert_with_parent(collider, handle, &mut bodies); let bo_handle = fluids_pipeline .liquid_world @@ -108,7 +110,7 @@ pub fn init_world(testbed: &mut Testbed) { let samples = salva2d::sampling::shape_surface_ray_sample(collider.shape(), PARTICLE_RADIUS).unwrap(); let rb = RigidBodyBuilder::dynamic() - .translation(Vector2::new(x, y)) + .translation(Vector2::new(x, y).into()) .build(); let rb_handle = bodies.insert(rb); let co_handle = colliders.insert_with_parent(collider, rb_handle, &mut bodies); @@ -133,15 +135,16 @@ pub fn init_world(testbed: &mut Testbed) { * Set up the testbed. */ plugin.set_pipeline(fluids_pipeline); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; + testbed.look_at(Vector2::new(0.0, 5.5).into(), 50.0); // testbed.enable_boundary_particles_rendering(true); } diff --git a/examples2d/custom_forces2.rs b/examples2d/custom_forces2.rs index a97e5ab..5bf8cf7 100644 --- a/examples2d/custom_forces2.rs +++ b/examples2d/custom_forces2.rs @@ -1,6 +1,6 @@ extern crate nalgebra as na; -use na::{Point2, Point3, Unit, Vector2}; +use na::{Unit, Vector2, Vector3}; use rapier2d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodySet}; use rapier2d::geometry::ColliderSet; use rapier_testbed2d::Testbed; @@ -30,37 +30,37 @@ pub fn init_world(testbed: &mut Testbed) { // Liquid. let nparticles = 30; let custom_force1 = CustomForceField { - origin: Point2::new(1.0, 0.0), + origin: Vector2::new(1.0, 0.0), }; let custom_force2 = CustomForceField { - origin: Point2::new(-1.0, 0.0), + origin: Vector2::new(-1.0, 0.0), }; let mut fluid = helper::cube_fluid(nparticles, nparticles, PARTICLE_RADIUS, 1000.0); fluid.nonpressure_forces.push(Box::new(custom_force1)); fluid.nonpressure_forces.push(Box::new(custom_force2)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); /* * Set up the testbed. */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point2::origin(), 300.0); + testbed.look_at(Vector2::zeros().into(), 300.0); } struct CustomForceField { - origin: Point2, + origin: Vector2, } impl NonPressureForce for CustomForceField { diff --git a/examples2d/elasticity2.rs b/examples2d/elasticity2.rs index 147df00..67dc2d8 100644 --- a/examples2d/elasticity2.rs +++ b/examples2d/elasticity2.rs @@ -1,6 +1,6 @@ extern crate nalgebra as na; -use na::{Isometry2, Point2, Point3, Vector2}; +use na::{Isometry2, Vector2, Vector3}; use rapier2d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier2d::geometry::{ColliderBuilder, ColliderSet}; use rapier_testbed2d::Testbed; @@ -49,7 +49,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); // Second fluid with smaller young modulus. let elasticity: Becker2009Elasticity = Becker2009Elasticity::new(100_000.0, 0.3, true); @@ -61,7 +61,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.6, 0.8, 0.5)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.6, 0.8, 0.5)); // Setup the ground. let ground_handle = bodies.insert(RigidBodyBuilder::fixed().build()); @@ -81,15 +81,15 @@ pub fn init_world(testbed: &mut Testbed) { */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point2::new(0.0, 1.0), 100.0); + testbed.look_at(Vector2::new(0.0, 1.0).into(), 100.0); } diff --git a/examples2d/helper.rs b/examples2d/helper.rs index 192c1ed..25a6137 100644 --- a/examples2d/helper.rs +++ b/examples2d/helper.rs @@ -1,4 +1,4 @@ -use na::{Point2, Vector2}; +use na::Vector2; use salva2d::object::{interaction_groups::InteractionGroups, Fluid}; pub fn cube_fluid(ni: usize, nj: usize, particle_rad: f32, density: f32) -> Fluid { @@ -9,7 +9,7 @@ pub fn cube_fluid(ni: usize, nj: usize, particle_rad: f32, density: f32) -> Flui for j in 0..nj { let x = (i as f32) * particle_rad * 2.0; let y = (j as f32) * particle_rad * 2.0; - points.push(Point2::new(x, y) + Vector2::repeat(particle_rad) - half_extents); + points.push(Vector2::new(x, y) + Vector2::repeat(particle_rad) - half_extents); } } diff --git a/examples2d/layers2.rs b/examples2d/layers2.rs index 20155aa..d92c4f1 100644 --- a/examples2d/layers2.rs +++ b/examples2d/layers2.rs @@ -1,8 +1,8 @@ extern crate nalgebra as na; -use na::{DVector, Point2, Point3, Vector2}; +use na::{Vector2, Vector3}; use rapier2d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; -use rapier2d::geometry::{Collider, ColliderBuilder, ColliderSet}; +use rapier2d::geometry::{Collider, ColliderBuilder, ColliderSet, InteractionTestMode}; use rapier_testbed2d::Testbed; use salva2d::integrations::rapier::{ColliderSampling, FluidsPipeline, FluidsTestbedPlugin}; use salva2d::object::interaction_groups::{Group, InteractionGroups}; @@ -38,8 +38,8 @@ pub fn init_world(testbed: &mut Testbed) { for j in 0..nj { let x = (i as f32) * PARTICLE_RADIUS * 2.0 - ni as f32 * PARTICLE_RADIUS; let y = (j as f32 + 1.0) * PARTICLE_RADIUS * 2.0 + 0.5; - points1.push(Point2::new(x, y)); - points2.push(Point2::new(x + ni as f32 * PARTICLE_RADIUS, y)); + points1.push(Vector2::new(x, y)); + points2.push(Vector2::new(x + ni as f32 * PARTICLE_RADIUS, y)); } } @@ -47,7 +47,7 @@ pub fn init_world(testbed: &mut Testbed) { for j in 0..nj * 2 { let x = (i as f32) * PARTICLE_RADIUS * 2.0 - ni as f32 * PARTICLE_RADIUS; let y = (j as f32 + 1.0) * PARTICLE_RADIUS * 2.0 + 0.5; - points3.push(Point2::new(x, y + shift2)); + points3.push(Vector2::new(x, y + shift2)); } } @@ -62,7 +62,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); let elasticity: Becker2009Elasticity = Becker2009Elasticity::new(1_000.0, 0.3, true); let viscosity = XSPHViscosity::new(0.5, 1.0); @@ -75,7 +75,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(1.0, 0.4, 0.6)); + plugin.set_fluid_color(fluid_handle, Vector3::new(1.0, 0.4, 0.6)); let viscosity = ArtificialViscosity::new(0.5, 0.0); let mut fluid = Fluid::new( @@ -86,7 +86,7 @@ pub fn init_world(testbed: &mut Testbed) { ); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.6, 0.8, 0.5)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.6, 0.8, 0.5)); /* * Ground @@ -94,17 +94,19 @@ pub fn init_world(testbed: &mut Testbed) { let ground_size = Vector2::new(10.0, 1.0); let nsubdivs = 50; - let heights = DVector::from_fn(nsubdivs + 1, |i, _| { - if i == 0 || i == nsubdivs { - 20.0 - } else { - (i as f32 * ground_size.x / (nsubdivs as f32)).cos() * 0.5 - } - }); + let heights: Vec<_> = (0..=nsubdivs) + .map(|i| { + if i == 0 || i == nsubdivs { + 20.0 + } else { + (i as f32 * ground_size.x / (nsubdivs as f32)).cos() * 0.5 + } + }) + .collect(); let rigid_body = RigidBodyBuilder::fixed().build(); let handle = bodies.insert(rigid_body); - let collider = ColliderBuilder::heightfield(heights, ground_size).build(); + let collider = ColliderBuilder::heightfield(heights, ground_size.into()).build(); let co_handle = colliders.insert_with_parent(collider, handle, &mut bodies); let bo_handle = fluids_pipeline .liquid_world @@ -125,7 +127,7 @@ pub fn init_world(testbed: &mut Testbed) { salva2d::sampling::shape_surface_ray_sample(collider.shape(), PARTICLE_RADIUS) .unwrap(); let rb = RigidBodyBuilder::dynamic() - .translation(Vector2::new(x, y)) + .translation(Vector2::new(x, y).into()) .build(); let rb_handle = bodies.insert(rb); let membership: u32 = interaction_group.memberships.into(); @@ -133,6 +135,7 @@ pub fn init_world(testbed: &mut Testbed) { collider.set_collision_groups(rapier2d::geometry::InteractionGroups::new( rapier2d::geometry::Group::from(membership), rapier2d::geometry::Group::from(filter), + InteractionTestMode::And, )); let co_handle = colliders.insert_with_parent(collider, rb_handle, &mut bodies); let bo_handle = fluids_pipeline @@ -171,15 +174,16 @@ pub fn init_world(testbed: &mut Testbed) { * Set up the testbed. */ plugin.set_pipeline(fluids_pipeline); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; + testbed.look_at(Vector2::new(0.0, 5.5).into(), 50.0); // testbed.enable_boundary_particles_rendering(true); } diff --git a/examples2d/surface_tension2.rs b/examples2d/surface_tension2.rs index 3d2549c..d677a4d 100644 --- a/examples2d/surface_tension2.rs +++ b/examples2d/surface_tension2.rs @@ -1,9 +1,9 @@ extern crate nalgebra as na; -use na::{Isometry2, Point2, Point3, Vector2}; +use na::{Isometry2, Vector2, Vector3}; use rapier2d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier2d::geometry::{ColliderBuilder, ColliderSet}; -use rapier_testbed2d::{Testbed, TestbedApp}; +use rapier_testbed2d::{Example, Testbed, TestbedApp}; use salva2d::integrations::rapier::{ ColliderSampling, FluidsPipeline, FluidsRenderingMode, FluidsTestbedPlugin, }; @@ -40,7 +40,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(surface_tension)); fluid.nonpressure_forces.push(Box::new(viscosity)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); // Setup the ground. let ground_thickness = 0.02; @@ -63,21 +63,21 @@ pub fn init_world(testbed: &mut Testbed) { */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point2::origin(), 1500.0); + testbed.look_at(Vector2::zeros().into(), 1500.0); // testbed.enable_boundary_particles_rendering(true); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Boxes", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Boxes", init_world)]); + pollster::block_on(testbed.run()); } diff --git a/examples3d/Cargo.toml b/examples3d/Cargo.toml index 650c9cf..5d84140 100644 --- a/examples3d/Cargo.toml +++ b/examples3d/Cargo.toml @@ -10,12 +10,13 @@ parallel = ["rapier_testbed3d/parallel", "salva3d/parallel"] [dependencies] num-traits = "0.2" -Inflector = "0.11" -nalgebra = "0.33" -rapier3d = "0.23" -rapier_testbed3d = "0.23" -parry3d = "0.18" +Inflector = "0.11" +nalgebra = "0" +rapier3d = "0.33" +rapier_testbed3d = "0.33" +parry3d = "0.28" bevy = "0.15" +pollster = "0.4" [dependencies.salva3d] path = "../build/salva3d" diff --git a/examples3d/all_examples3.rs b/examples3d/all_examples3.rs index 15afac7..bd28e41 100644 --- a/examples3d/all_examples3.rs +++ b/examples3d/all_examples3.rs @@ -5,7 +5,7 @@ use wasm_bindgen::prelude::*; use inflector::Inflector; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, TestbedApp}; use std::cmp::Ordering; mod basic3; @@ -51,27 +51,30 @@ pub fn main() { .unwrap_or(String::new()) .to_camel_case(); - let mut builders: Vec<(_, fn(&mut Testbed))> = vec![ - ("Basic", basic3::init_world), - ("Height field", heightfield3::init_world), - ("Custom Forces", custom_forces3::init_world), - ("Elasticity", elasticity3::init_world), - ("Faucet", faucet3::init_world), //FIXME: bug with adding & removing particles - ("Surface tension", surface_tension3::init_world), + let mut builders = vec![ + Example::demo("Basic", basic3::init_world), + Example::demo("Height field", heightfield3::init_world), + Example::demo("Custom Forces", custom_forces3::init_world), + Example::demo("Elasticity", elasticity3::init_world), + Example::demo("Faucet", faucet3::init_world), //FIXME: bug with adding & removing particles + Example::demo("Surface tension", surface_tension3::init_world), ]; // Lexicographic sort, with stress tests moved at the end of the list. - builders.sort_by(|a, b| match (a.0.starts_with("("), b.0.starts_with("(")) { - (true, true) | (false, false) => a.0.cmp(b.0), - (true, false) => Ordering::Greater, - (false, true) => Ordering::Less, - }); + builders.sort_by( + |a, b| match (a.name.starts_with("("), b.name.starts_with("(")) { + (true, true) | (false, false) => a.name.cmp(b.name), + (true, false) => Ordering::Greater, + (false, true) => Ordering::Less, + }, + ); let i = builders .iter() - .position(|builder| builder.0.to_camel_case().as_str() == demo.as_str()) + .position(|builder| builder.name.to_camel_case().as_str() == demo.as_str()) .unwrap_or(0); + builders.rotate_left(i); - let testbed = TestbedApp::from_builders(i, builders); - testbed.run() + let testbed = TestbedApp::from_builders(builders); + pollster::block_on(testbed.run()); } diff --git a/examples3d/basic3.rs b/examples3d/basic3.rs index 5423195..a99dcc2 100644 --- a/examples3d/basic3.rs +++ b/examples3d/basic3.rs @@ -1,9 +1,9 @@ extern crate nalgebra as na; -use na::{Isometry3, Point3, Vector3}; +use na::{Isometry3, Vector3}; use rapier3d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier3d::geometry::{ColliderBuilder, ColliderSet, SharedShape}; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, Testbed, TestbedApp}; use salva3d::integrations::rapier::{ColliderSampling, FluidsPipeline, FluidsTestbedPlugin}; use salva3d::object::interaction_groups::InteractionGroups; use salva3d::object::Boundary; @@ -70,7 +70,7 @@ pub fn init_world(testbed: &mut Testbed) { let samples = salva3d::sampling::shape_surface_ray_sample(&*wall_shape, PARTICLE_RADIUS).unwrap(); let co = ColliderBuilder::new(wall_shape.clone()) - .position(*pose) + .position((*pose).into()) .build(); let co_handle = colliders.insert_with_parent(co, ground_handle, &mut bodies); let bo_handle = fluids_pipeline @@ -103,23 +103,23 @@ pub fn init_world(testbed: &mut Testbed) { */ let mut plugin = FluidsTestbedPlugin::new(); plugin.set_pipeline(fluids_pipeline); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); plugin.render_boundary_particles = true; - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); // testbed.set_body_wireframe(ground_handle, true); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point3::new(3.0, 3.0, 3.0), Point3::origin()); + testbed.look_at(Vector3::new(3.0, 3.0, 3.0).into(), Vector3::zeros().into()); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Basic", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Basic", init_world)]); + pollster::block_on(testbed.run()); } diff --git a/examples3d/custom_forces3.rs b/examples3d/custom_forces3.rs index 2158d38..7d4f4f3 100644 --- a/examples3d/custom_forces3.rs +++ b/examples3d/custom_forces3.rs @@ -1,9 +1,9 @@ extern crate nalgebra as na; -use na::{Point3, Unit, Vector3}; +use na::{Unit, Vector3}; use rapier3d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodySet}; use rapier3d::geometry::ColliderSet; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, Testbed, TestbedApp}; use salva3d::integrations::rapier::{FluidsPipeline, FluidsRenderingMode, FluidsTestbedPlugin}; use salva3d::object::{Boundary, Fluid}; use salva3d::solver::NonPressureForce; @@ -30,42 +30,42 @@ pub fn init_world(testbed: &mut Testbed) { // fluids. let nparticles = 10; let custom_force1 = CustomForceField { - origin: Point3::new(1.0, 0.0, 0.0), + origin: Vector3::new(1.0, 0.0, 0.0), }; let custom_force2 = CustomForceField { - origin: Point3::new(-1.0, 0.0, 0.0), + origin: Vector3::new(-1.0, 0.0, 0.0), }; let mut fluid = helper::cube_fluid(nparticles, nparticles, nparticles, PARTICLE_RADIUS, 1000.0); fluid.nonpressure_forces.push(Box::new(custom_force1)); fluid.nonpressure_forces.push(Box::new(custom_force2)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); /* * Set up the testbed. */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point3::new(3.0, 3.0, 3.0), Point3::origin()); + testbed.look_at(Vector3::new(3.0, 3.0, 3.0).into(), Vector3::zeros().into()); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Boxes", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Boxes", init_world)]); + pollster::block_on(testbed.run()); } struct CustomForceField { - origin: Point3, + origin: Vector3, } impl NonPressureForce for CustomForceField { diff --git a/examples3d/elasticity3.rs b/examples3d/elasticity3.rs index 1ad8198..aefa1d1 100644 --- a/examples3d/elasticity3.rs +++ b/examples3d/elasticity3.rs @@ -1,9 +1,9 @@ extern crate nalgebra as na; -use na::{Isometry3, Point3, Vector3}; +use na::{Isometry3, Vector3}; use rapier3d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier3d::geometry::{ColliderBuilder, ColliderSet}; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, Testbed, TestbedApp}; use salva3d::integrations::rapier::{ ColliderSampling, FluidsPipeline, FluidsRenderingMode, FluidsTestbedPlugin, }; @@ -56,7 +56,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity.clone())); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); // Second fluid with smaller young modulus. let elasticity: Becker2009Elasticity = Becker2009Elasticity::new(100_000.0, 0.3, true); @@ -75,7 +75,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(elasticity)); fluid.nonpressure_forces.push(Box::new(viscosity)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.6, 0.8, 0.5)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.6, 0.8, 0.5)); // Setup the ground. let ground_handle = bodies.insert(RigidBodyBuilder::fixed().build()); @@ -96,21 +96,21 @@ pub fn init_world(testbed: &mut Testbed) { */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_body_wireframe(ground_handle, true); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point3::new(1.5, 1.5, 1.5), Point3::origin()); + testbed.look_at(Vector3::new(1.5, 1.5, 1.5).into(), Vector3::zeros().into()); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Elasticity", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Elasticity", init_world)]); + pollster::block_on(testbed.run()); } diff --git a/examples3d/faucet3.rs b/examples3d/faucet3.rs index dc2e5ec..d20912f 100644 --- a/examples3d/faucet3.rs +++ b/examples3d/faucet3.rs @@ -1,12 +1,12 @@ extern crate nalgebra as na; -use na::{Point3, Vector3}; +use na::Vector3; use rapier3d::geometry::{ColliderBuilder, ColliderSet}; use rapier3d::{ dynamics::{RigidBodyBuilder, RigidBodySet}, prelude::{ImpulseJointSet, MultibodyJointSet}, }; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, Testbed, TestbedApp}; use salva3d::integrations::rapier::{ColliderSampling, FluidsPipeline, FluidsTestbedPlugin}; use salva3d::object::interaction_groups::InteractionGroups; use salva3d::object::{Boundary, Fluid}; @@ -45,7 +45,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(viscosity)); fluid.nonpressure_forces.push(Box::new(tension)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.5, 1.0, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.5, 1.0, 1.0)); // Setup the ground. let ground_handle = bodies.insert(RigidBodyBuilder::fixed().build()); @@ -95,7 +95,7 @@ pub fn init_world(testbed: &mut Testbed) { for i in 0..nparticles { for j in 0..nparticles { - let pos = Point3::new(i as f32 * diam, height, j as f32 * diam); + let pos = Vector3::new(i as f32 * diam, height, j as f32 * diam); particles.push(pos + Vector3::new(shift, 0.0, shift)); velocities.push(Vector3::y() * vel); } @@ -108,22 +108,25 @@ pub fn init_world(testbed: &mut Testbed) { * Set up the testbed. */ plugin.set_pipeline(fluids_pipeline); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_body_wireframe(ground_handle, true); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; // testbed.enable_boundary_particles_rendering(true); - testbed.look_at(Point3::new(1.5, 0.0, 1.5), Point3::new(0.0, 0.0, 0.0)); + testbed.look_at( + Vector3::new(1.5, 0.0, 1.5).into(), + Vector3::new(0.0, 0.0, 0.0).into(), + ); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Boxes", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Boxes", init_world)]); + pollster::block_on(testbed.run()); } diff --git a/examples3d/harness_basic3.rs b/examples3d/harness_basic3.rs index 65d9abc..7edc572 100644 --- a/examples3d/harness_basic3.rs +++ b/examples3d/harness_basic3.rs @@ -6,7 +6,7 @@ use rapier3d::{ dynamics::{RigidBodyBuilder, RigidBodySet}, prelude::{ImpulseJointSet, MultibodyJointSet}, }; -use rapier_testbed3d::harness::Harness; +use rapier_testbed3d::harness::{Harness, RapierBroadPhaseType}; use salva3d::integrations::rapier::{ColliderSampling, FluidsHarnessPlugin, FluidsPipeline}; use salva3d::object::interaction_groups::InteractionGroups; use salva3d::object::Boundary; @@ -74,7 +74,7 @@ pub fn init_world(harness: &mut Harness) { let samples = salva3d::sampling::shape_surface_ray_sample(&*wall_shape, PARTICLE_RADIUS).unwrap(); let co = ColliderBuilder::new(wall_shape.clone()) - .position(*pose) + .position((*pose).into()) .build(); let co_handle = colliders.insert_with_parent(co, ground_handle, &mut bodies); let bo_handle = fluids_pipeline @@ -113,7 +113,8 @@ pub fn init_world(harness: &mut Harness) { colliders, impulse_joints, multibody_joints, - gravity, + RapierBroadPhaseType::default(), + gravity.into(), (), ); harness.integration_parameters_mut().dt = 1.0 / 200.0; diff --git a/examples3d/heightfield3.rs b/examples3d/heightfield3.rs index e09b56c..32d8cdd 100644 --- a/examples3d/heightfield3.rs +++ b/examples3d/heightfield3.rs @@ -1,6 +1,7 @@ extern crate nalgebra as na; -use nalgebra::Isometry3; +use nalgebra::{Isometry3, Vector3}; +use rapier3d::geometry::Array2; use rapier3d::na::ComplexField; use rapier3d::prelude::*; use rapier_testbed3d::Testbed; @@ -37,7 +38,7 @@ pub fn init_world(testbed: &mut Testbed) { )); let viscosity = ArtificialViscosity::new(1.0, 0.0); fluid.nonpressure_forces.push(Box::new(viscosity)); - fluid.velocities = vec![-Vector::y() * 10.; fluid.velocities.len()]; + fluid.velocities = vec![-Vector3::y() * 10.; fluid.velocities.len()]; let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); /* @@ -46,19 +47,27 @@ pub fn init_world(testbed: &mut Testbed) { let ground_size = Vector::new(12.0, 1.0, 12.0); let nsubdivs = 40; - let heights = DMatrix::from_fn(nsubdivs + 1, nsubdivs + 1, |i, j| { - if i == 0 || i == nsubdivs || j == 0 || j == nsubdivs { - 3.0 - } else { - let x = i as f32 * ground_size.x / (nsubdivs as f32); - let z = j as f32 * ground_size.z / (nsubdivs as f32); + let heights = Array2::new( + nsubdivs + 1, + nsubdivs + 1, + (0..=nsubdivs) + .flat_map(|j| { + (0..=nsubdivs).map(move |i| { + if i == 0 || i == nsubdivs || j == 0 || j == nsubdivs { + 3.0 + } else { + let x = i as f32 * ground_size.x / (nsubdivs as f32); + let z = j as f32 * ground_size.z / (nsubdivs as f32); - // NOTE: make sure we use the sin/cos from simba to ensure - // cross-platform determinism of the example when the - // enhanced_determinism feature is enabled. - ComplexField::sin(x) + ComplexField::cos(z) - } - }); + // NOTE: make sure we use the sin/cos from simba to ensure + // cross-platform determinism of the example when the + // enhanced_determinism feature is enabled. + ComplexField::sin(x) + ComplexField::cos(z) + } + }) + }) + .collect(), + ); let rigid_body = RigidBodyBuilder::fixed().build(); let handle = bodies.insert(rigid_body); @@ -81,15 +90,15 @@ pub fn init_world(testbed: &mut Testbed) { let mut plugin = FluidsTestbedPlugin::new(); plugin.set_pipeline(fluids_pipeline); - plugin.set_fluid_color(fluid_handle, Point::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector::new(0.8, 0.7, 1.0).into()); // plugin.render_boundary_particles = true; - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_body_wireframe(handle, true); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - // testbed.look_at(Point3::new(3.0, 3.0, 3.0), Point3::origin()); + // testbed.look_at(Vector3::new(3.0, 3.0, 3.0), Vector3::origin()); /* * Set up the testbed. */ testbed.set_world(bodies, colliders, impulse_joints, multibody_joints); - testbed.look_at(point![100.0, 100.0, 100.0], Point::origin()); + testbed.look_at(Vector::new(100.0, 100.0, 100.0), Vector::ZERO); } diff --git a/examples3d/helper.rs b/examples3d/helper.rs index 5735661..8426060 100644 --- a/examples3d/helper.rs +++ b/examples3d/helper.rs @@ -1,4 +1,4 @@ -use super::na::{Point3, Vector3}; +use super::na::Vector3; use salva3d::object::{interaction_groups::InteractionGroups, Fluid}; pub fn cube_fluid(ni: usize, nj: usize, nk: usize, particle_rad: f32, density: f32) -> Fluid { @@ -11,7 +11,7 @@ pub fn cube_fluid(ni: usize, nj: usize, nk: usize, particle_rad: f32, density: f let x = (i as f32) * particle_rad * 2.0; let y = (j as f32) * particle_rad * 2.0; let z = (k as f32) * particle_rad * 2.0; - points.push(Point3::new(x, y, z) + Vector3::repeat(particle_rad) - half_extents); + points.push(Vector3::new(x, y, z) + Vector3::repeat(particle_rad) - half_extents); } } } diff --git a/examples3d/surface_tension3.rs b/examples3d/surface_tension3.rs index 7ad5413..292aa9f 100644 --- a/examples3d/surface_tension3.rs +++ b/examples3d/surface_tension3.rs @@ -1,9 +1,9 @@ extern crate nalgebra as na; -use na::{Isometry3, Point3, Vector3}; +use na::{Isometry3, Vector3}; use rapier3d::dynamics::{ImpulseJointSet, MultibodyJointSet, RigidBodyBuilder, RigidBodySet}; use rapier3d::geometry::{ColliderBuilder, ColliderSet}; -use rapier_testbed3d::{Testbed, TestbedApp}; +use rapier_testbed3d::{Example, Testbed, TestbedApp}; use salva3d::integrations::rapier::{ ColliderSampling, FluidsPipeline, FluidsRenderingMode, FluidsTestbedPlugin, }; @@ -43,7 +43,7 @@ pub fn init_world(testbed: &mut Testbed) { fluid.nonpressure_forces.push(Box::new(surface_tension)); fluid.nonpressure_forces.push(Box::new(viscosity)); let fluid_handle = fluids_pipeline.liquid_world.add_fluid(fluid); - plugin.set_fluid_color(fluid_handle, Point3::new(0.8, 0.7, 1.0)); + plugin.set_fluid_color(fluid_handle, Vector3::new(0.8, 0.7, 1.0)); // Setup the ground. let ground_handle = bodies.insert(RigidBodyBuilder::fixed().build()); @@ -69,21 +69,24 @@ pub fn init_world(testbed: &mut Testbed) { */ plugin.set_pipeline(fluids_pipeline); plugin.set_fluid_rendering_mode(FluidsRenderingMode::VelocityColor { min: 0.0, max: 5.0 }); - testbed.add_plugin(plugin); + plugin.add_to_testbed(testbed); testbed.set_body_wireframe(ground_handle, true); testbed.set_world_with_params( bodies, colliders, impulse_joints, multibody_joints, - gravity, + gravity.into(), (), ); testbed.integration_parameters_mut().dt = 1.0 / 200.0; - testbed.look_at(Point3::new(0.25, 0.25, 0.25), Point3::origin()); + testbed.look_at( + Vector3::new(0.25, 0.25, 0.25).into(), + Vector3::zeros().into(), + ); } fn main() { - let testbed = TestbedApp::from_builders(0, vec![("Surface tension", init_world)]); - testbed.run() + let testbed = TestbedApp::from_builders(vec![Example::demo("Surface tension", init_world)]); + pollster::block_on(testbed.run()); } diff --git a/src/counters/timer.rs b/src/counters/timer.rs index 9add9d3..b2a03be 100644 --- a/src/counters/timer.rs +++ b/src/counters/timer.rs @@ -37,15 +37,15 @@ impl Timer { pub fn start(&mut self) { if self.enabled { self.time = 0.0; - self.start = Some(instant::now()); + //self.start = Some(instant::now()); } } /// Pause the timer. pub fn pause(&mut self) { if self.enabled { - if let Some(start) = self.start { - self.time += instant::now() - start; + if self.start.is_some() { + //self.time += instant::now() - start; } self.start = None; } @@ -54,7 +54,7 @@ impl Timer { /// Resume the timer. pub fn resume(&mut self) { if self.enabled { - self.start = Some(instant::now()); + //self.start = Some(instant::now()); } } diff --git a/src/coupling/coupling_manager.rs b/src/coupling/coupling_manager.rs index 3eeea2c..20e9a40 100644 --- a/src/coupling/coupling_manager.rs +++ b/src/coupling/coupling_manager.rs @@ -24,7 +24,12 @@ pub trait CouplingManager { ); /// Transmit forces from salva's boundary objects to the coupled bodies. - fn transmit_forces(&mut self, timestep: &TimestepManager, boundaries: &BoundarySet); + fn transmit_forces( + &mut self, + timestep: &TimestepManager, + boundaries: &BoundarySet, + force_coefficient: Real, + ); } impl CouplingManager for () { @@ -39,5 +44,5 @@ impl CouplingManager for () { ) { } - fn transmit_forces(&mut self, _: &TimestepManager, _: &BoundarySet) {} + fn transmit_forces(&mut self, _: &TimestepManager, _: &BoundarySet, _: Real) {} } diff --git a/src/geometry/contacts.rs b/src/geometry/contacts.rs index a4d0964..6cf2f04 100644 --- a/src/geometry/contacts.rs +++ b/src/geometry/contacts.rs @@ -1,6 +1,6 @@ use crate::counters::Counters; use crate::geometry::HGrid; -use crate::math::{Point, Real, Vector}; +use crate::math::{Real, Vector}; use crate::object::Boundary; use crate::object::Fluid; @@ -258,9 +258,9 @@ fn compute_contacts_for_pair_of_cells( fluid_fluid_contacts: &[ParticlesContacts], fluid_boundary_contacts: &[ParticlesContacts], boundary_boundary_contacts: &[ParticlesContacts], - curr_cell: &Point, + curr_cell: &Vector, curr_particles: &[HGridEntry], - neighbor_cell: &Point, + neighbor_cell: &Vector, neighbor_particles: &[HGridEntry], ) { for entry in curr_particles { @@ -282,7 +282,7 @@ fn compute_contacts_for_pair_of_cells( let pi = &bi.positions[*particle_i]; let pj = &bj.positions[*particle_j]; - if na::distance_squared(pi, pj) <= h * h { + if (*pi - *pj).norm_squared() <= h * h { let contact = Contact { i_model: *boundary_i, j_model: *boundary_j, @@ -319,7 +319,7 @@ fn compute_contacts_for_pair_of_cells( let pi = &boundaries[*boundary_i].positions[*particle_i]; let pj = &fluids[*fluid_j].positions[*particle_j]; - if na::distance_squared(pi, pj) <= h * h { + if (*pi - *pj).norm_squared() <= h * h { let contact = Contact { i_model: *fluid_j, j_model: *boundary_i, @@ -363,8 +363,8 @@ fn compute_contacts_for_pair_of_cells( fluids[fluid_j].positions[particle_j] }; - if na::distance_squared(&pi, &pj) <= h * h { - assert!(na::distance_squared(&pj, &pi) <= h * h); + if (pi - pj).norm_squared() <= h * h { + assert!((pj - pi).norm_squared() <= h * h); let contact = Contact { i_model: *fluid_i, j_model: fluid_j, @@ -424,7 +424,7 @@ pub fn compute_self_contacts(h: Real, fluid: &Fluid, contacts: &mut ParticlesCon let pi = fluid.positions[*particle_i]; let pj = fluid.positions[*particle_j]; - if na::distance_squared(&pi, &pj) <= h * h { + if (pi - pj).norm_squared() <= h * h { let contact = Contact { i_model: 0, j_model: 0, diff --git a/src/geometry/hgrid.rs b/src/geometry/hgrid.rs index 4abab2a..531434d 100644 --- a/src/geometry/hgrid.rs +++ b/src/geometry/hgrid.rs @@ -1,8 +1,7 @@ use fnv::FnvHasher; - use std::collections::HashMap; -use crate::math::{Point, Real, Vector, DIM}; +use crate::math::{Real, Vector, DIM}; use std::hash::BuildHasher; @@ -20,7 +19,7 @@ impl BuildHasher for DeterministicState { /// A grid based on spacial hashing. #[derive(PartialEq, Debug, Clone)] pub struct HGrid { - cells: HashMap, Vec, DeterministicState>, + cells: HashMap, Vec, DeterministicState>, cell_width: Real, } @@ -47,8 +46,8 @@ impl HGrid { } /// Computes the logical grid cell containing `point`. - pub fn key(&self, point: &Point) -> Point { - Point::from(point.coords.map(|e| Self::quantify(e, self.cell_width))) + pub fn key(&self, point: &Vector) -> Vector { + Vector::from(point.map(|e| Self::quantify(e, self.cell_width))) } /// Removes all elements from this grid. @@ -57,7 +56,7 @@ impl HGrid { } /// Inserts the given `element` into the cell containing the given `point`. - pub fn insert(&mut self, point: &Point, element: T) { + pub fn insert(&mut self, point: &Vector, element: T) { let key = self.key(point); self.cells.entry(key).or_insert(Vec::new()).push(element) } @@ -65,7 +64,7 @@ impl HGrid { /// Returns the element attached to the cell containing the given `point`. /// /// Returns `None` if the cell is empty. - pub fn cell_containing_point(&self, point: &Point) -> Option<&Vec> { + pub fn cell_containing_point(&self, point: &Vector) -> Option<&Vec> { let key = self.key(point); self.cells.get(&key) } @@ -73,17 +72,17 @@ impl HGrid { /// An iterator through all the non-empty cells of this grid. /// /// The returned tuple include the cell indentifier, and the elements attached to this cell. - pub fn cells(&self) -> impl Iterator, &Vec)> { + pub fn cells(&self) -> impl Iterator, &Vec)> { self.cells.iter() } /// The underlying hash map of this spacial grid. - pub fn inner_table(&self) -> &HashMap, Vec, DeterministicState> { + pub fn inner_table(&self) -> &HashMap, Vec, DeterministicState> { &self.cells } /// Get the content of the logical cell identified by `key`. - pub fn cell(&self, key: &Point) -> Option<&Vec> { + pub fn cell(&self, key: &Vector) -> Option<&Vec> { self.cells.get(key) } @@ -92,9 +91,9 @@ impl HGrid { /// The given cell itself will be yielded by this iterator too. pub fn neighbor_cells( &self, - cell: &Point, + cell: &Vector, radius: Real, - ) -> impl Iterator, &Vec)> { + ) -> impl Iterator, &Vec)> { let cells = &self.cells; let quantified_radius = Self::quantify_ceil(radius, self.cell_width); @@ -104,7 +103,7 @@ impl HGrid { // pub fn elements_close_to_point<'a>( // &'a self, - // point: &Point, + // point: &Vector, // radius: Real, // ) -> impl Iterator // { @@ -121,9 +120,9 @@ impl HGrid { /// An iterator through all the cells intersecting the given AABB. pub fn cells_intersecting_aabb( &self, - mins: &Point, - maxs: &Point, - ) -> impl Iterator, &Vec)> { + mins: &Vector, + maxs: &Vector, + ) -> impl Iterator, &Vec)> { let cells = &self.cells; let start = self.key(mins); let end = self.key(maxs); @@ -132,20 +131,20 @@ impl HGrid { .filter_map(move |cell| cells.get(&cell).map(|c| (cell, c))) } - // pub fn elements_containing_point(&self, point: &Point) -> impl Iterator { + // pub fn elements_containing_point(&self, point: &Vector) -> impl Iterator { // std::iter::empty() // } } struct CellRangeIterator { - start: Point, - end: Point, - curr: Point, + start: Vector, + end: Vector, + curr: Vector, done: bool, } impl CellRangeIterator { - fn new(start: Point, end: Point) -> Self { + fn new(start: Vector, end: Vector) -> Self { Self { start, end, @@ -154,7 +153,7 @@ impl CellRangeIterator { } } - fn with_center(center: Point, radius: i64) -> Self { + fn with_center(center: Vector, radius: i64) -> Self { let start = center - Vector::repeat(radius as i64); Self { start, @@ -166,7 +165,7 @@ impl CellRangeIterator { } impl Iterator for CellRangeIterator { - type Item = Point; + type Item = Vector; fn next(&mut self) -> Option { if self.done { @@ -200,32 +199,32 @@ mod test { #[cfg(feature = "dim2")] fn grid_neighbor_iterator() { use super::CellRangeIterator; - use crate::math::Point; + use crate::math::Vector; let expected = [ - Point::new(-1, 0), - Point::new(0, 0), - Point::new(1, 0), - Point::new(2, 0), - Point::new(3, 0), - Point::new(-1, 1), - Point::new(0, 1), - Point::new(1, 1), - Point::new(2, 1), - Point::new(3, 1), - Point::new(-1, 2), - Point::new(0, 2), - Point::new(1, 2), - Point::new(2, 2), - Point::new(3, 2), - Point::new(-1, 3), - Point::new(0, 3), - Point::new(1, 3), - Point::new(2, 3), - Point::new(3, 3), + Vector::new(-1, 0), + Vector::new(0, 0), + Vector::new(1, 0), + Vector::new(2, 0), + Vector::new(3, 0), + Vector::new(-1, 1), + Vector::new(0, 1), + Vector::new(1, 1), + Vector::new(2, 1), + Vector::new(3, 1), + Vector::new(-1, 2), + Vector::new(0, 2), + Vector::new(1, 2), + Vector::new(2, 2), + Vector::new(3, 2), + Vector::new(-1, 3), + Vector::new(0, 3), + Vector::new(1, 3), + Vector::new(2, 3), + Vector::new(3, 3), ]; - let iter = CellRangeIterator::with_center(Point::new(1, 2), 2); + let iter = CellRangeIterator::with_center(Vector::new(1, 2), 2); assert!(iter.zip(expected.iter()).all(|(a, b)| a == *b)) } diff --git a/src/integrations/rapier/fluids_pipeline.rs b/src/integrations/rapier/fluids_pipeline.rs index 2406d4e..d5d5f9c 100644 --- a/src/integrations/rapier/fluids_pipeline.rs +++ b/src/integrations/rapier/fluids_pipeline.rs @@ -1,19 +1,39 @@ use crate::coupling::CouplingManager; use crate::geometry::{HGrid, HGridEntry}; use crate::object::{BoundaryHandle, BoundarySet, Fluid}; -use crate::solver::DFSPHSolver; -use crate::LiquidWorld; +use crate::solver::{DFSPHSolver, DfsphParameters, PressureSolver}; use crate::TimestepManager; +use crate::{math, LiquidWorld}; use approx::AbsDiffEq; use na::Unit; use rapier::dynamics::RigidBodySet; use rapier::geometry::{ColliderHandle, ColliderSet}; -use rapier::math::{Point, Vector}; +use rapier::math::Vector as RapierVector; use rapier::parry::bounding_volume::BoundingVolume; use rapier::parry::shape::FeatureId; use std::collections::HashMap; use std::sync::RwLock; +#[cfg(feature = "dim2")] +fn salva_to_rapier_vector(v: math::Vector) -> RapierVector { + RapierVector::new(v.x, v.y) +} + +#[cfg(feature = "dim3")] +fn salva_to_rapier_vector(v: math::Vector) -> RapierVector { + RapierVector::new(v.x, v.y, v.z) +} + +#[cfg(feature = "dim2")] +fn rapier_to_salva_vector(v: RapierVector) -> math::Vector { + math::Vector::new(v.x, v.y) +} + +#[cfg(feature = "dim3")] +fn rapier_to_salva_vector(v: RapierVector) -> math::Vector { + math::Vector::new(v.x, v.y, v.z) +} + /// Pipeline for particle-based fluid simulation. pub struct FluidsPipeline { // FIXME: keep these public? @@ -31,15 +51,85 @@ impl FluidsPipeline { /// - `particle_radius`: the radius of every particle for the fluid simulation. /// - `smoothing_factor`: the smoothing factor used to compute the SPH kernel radius. /// The kernel radius will be computed as `particle_radius * smoothing_factor * 2.0. - pub fn new(particle_radius: f32, smoothing_factor: f32) -> Self { + pub fn new(particle_radius: math::Real, smoothing_factor: math::Real) -> Self { + Self::new_with_boundary_coef(particle_radius, smoothing_factor, na::one::()) + } + + /// Initialize a new pipeline with custom DFSPH solver parameters. + pub fn new_with_dfsph_parameters( + particle_radius: math::Real, + smoothing_factor: math::Real, + dfsph_parameters: DfsphParameters, + ) -> Self { + let dfsph: DFSPHSolver = DFSPHSolver::with_parameters(dfsph_parameters); + Self::new_with_solver(dfsph, particle_radius, smoothing_factor) + } + + /// Initialize a new pipeline for fluids simulation with a custom boundary force coefficient. + /// + /// # Parameters + /// + /// - `particle_radius`: the radius of every particle for the fluid simulation. + /// - `smoothing_factor`: the smoothing factor used to compute the SPH kernel radius. + /// The kernel radius will be computed as `particle_radius * smoothing_factor * 2.0. + /// - `boundary_force_coefficient`: coefficient applied when transmitting forces from fluids to boundaries. + /// Use 1.0 for full force, 0.5 for half force, etc. + pub fn new_with_boundary_coef( + particle_radius: math::Real, + smoothing_factor: math::Real, + boundary_force_coefficient: math::Real, + ) -> Self { let dfsph: DFSPHSolver = DFSPHSolver::new(); + Self::new_with_solver_and_boundary_coef( + dfsph, + particle_radius, + smoothing_factor, + boundary_force_coefficient, + ) + } + /// Initialize a new pipeline with a custom pressure solver. + pub fn new_with_solver( + solver: impl PressureSolver + Send + Sync + 'static, + particle_radius: math::Real, + smoothing_factor: math::Real, + ) -> Self { + Self::new_with_solver_and_boundary_coef( + solver, + particle_radius, + smoothing_factor, + na::one::(), + ) + } + + /// Initialize a new pipeline with a custom pressure solver and boundary force coefficient. + pub fn new_with_solver_and_boundary_coef( + solver: impl PressureSolver + Send + Sync + 'static, + particle_radius: math::Real, + smoothing_factor: math::Real, + boundary_force_coefficient: math::Real, + ) -> Self { Self { - liquid_world: LiquidWorld::new(dfsph, particle_radius, smoothing_factor), + liquid_world: LiquidWorld::new( + solver, + particle_radius, + smoothing_factor, + boundary_force_coefficient, + ), coupling: ColliderCouplingSet::new(), } } + /// Current DFSPH tuning parameters, if this pipeline uses DFSPH. + pub fn dfsph_parameters(&self) -> Option { + self.liquid_world.dfsph_parameters() + } + + /// Set DFSPH tuning parameters if this pipeline uses DFSPH. + pub fn set_dfsph_parameters(&mut self, parameters: DfsphParameters) -> bool { + self.liquid_world.set_dfsph_parameters(parameters) + } + /// Advances the fluid simulation by `dt` seconds. /// /// All the fluid particles will be affected by an acceleration equal to `gravity`. @@ -47,14 +137,15 @@ impl FluidsPipeline { /// However, it will not integrate these forces. Use the `PhysicsPipeline` for this integration. pub fn step( &mut self, - gravity: &Vector, - dt: f32, + gravity: &RapierVector, + dt: math::Real, colliders: &ColliderSet, bodies: &mut RigidBodySet, ) { + let gravity = rapier_to_salva_vector(*gravity); self.liquid_world.step_with_coupling( dt, - gravity, + &gravity, &mut self.coupling.as_manager_mut(colliders, bodies), ) } @@ -66,7 +157,7 @@ pub enum ColliderSampling { /// /// It is recommended that those points are separated by a distance smaller or equal to twice /// the particle radius used to initialize the LiquidWorld. - StaticSampling(Vec>), + StaticSampling(Vec>), /// The collider shape is approximated by a dynamic set of points automatically computed based on contacts with fluid particles. DynamicContactSampling, } @@ -147,8 +238,8 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { fn update_boundaries( &mut self, timestep: &TimestepManager, - h: f32, - particle_radius: f32, + h: math::Real, + particle_radius: math::Real, hgrid: &HGrid, fluids: &mut [Fluid], boundaries: &mut BoundarySet, @@ -179,47 +270,66 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { match &coupling.sampling_method { ColliderSampling::StaticSampling(points) => { for pt in points { - boundary.positions.push(collider.position() * pt); - let velocity = body.map(|b| b.velocity_at_point(pt)); + let rapier_pt = salva_to_rapier_vector(*pt); + let world_pt = rapier_to_salva_vector(collider.position() * rapier_pt); + boundary.positions.push(world_pt); + let velocity = body.map(|b| b.velocity_at_point(rapier_pt)); - boundary - .velocities - .push(velocity.unwrap_or(Vector::zeros())); + boundary.velocities.push(rapier_to_salva_vector( + velocity.unwrap_or(RapierVector::ZERO), + )); } - boundary.volumes.resize(points.len(), na::zero::()); + boundary + .volumes + .resize(points.len(), na::zero::()); } ColliderSampling::DynamicContactSampling => { - let prediction = h * na::convert::<_, f32>(0.5); - let margin = particle_radius * na::convert::<_, f32>(0.1); + let prediction = h * na::convert::<_, math::Real>(0.5); + let margin = particle_radius * na::convert::<_, math::Real>(0.1); let collider_pos = collider.position(); let aabb = collider .shape() .compute_aabb(&collider_pos) .loosened(h + prediction); + let mins = rapier_to_salva_vector(aabb.mins); + let maxs = rapier_to_salva_vector(aabb.maxs); for particle in hgrid - .cells_intersecting_aabb(&aabb.mins, &aabb.maxs) + .cells_intersecting_aabb(&mins, &maxs) .flat_map(|e| e.1) { match particle { HGridEntry::FluidParticle(fluid_id, particle_id) => { let fluid = &mut fluids[*fluid_id]; + + // Check interaction groups between this fluid and the boundary. + // Use the fluid's interaction groups explicitly to mirror checks + // performed elsewhere in the codebase. + let fluid_groups = fluid.interaction_groups; + let boundary_groups = boundary.interaction_groups; + + if !fluid_groups.test(boundary_groups) { + continue; + } + let particle_pos = fluid.positions[*particle_id] + fluid.velocities[*particle_id] * timestep.dt(); + let rapier_particle_pos = salva_to_rapier_vector(particle_pos); - if aabb.contains_local_point(&particle_pos) { + if aabb.contains_local_point(rapier_particle_pos) { let (proj, feature) = collider.shape().project_point_and_get_feature( &collider_pos, - &particle_pos, + rapier_particle_pos, ); - let dpt = particle_pos - proj.point; + let dpt = particle_pos - rapier_to_salva_vector(proj.point); - if let Some((normal, depth)) = - Unit::try_new_and_get(dpt, f32::default_epsilon()) - { + if let Some((normal, depth)) = Unit::try_new_and_get( + dpt, + math::Real::default_epsilon(), + ) { if proj.is_inside { fluid.positions[*particle_id] -= *normal * (depth + margin); @@ -227,7 +337,7 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { let vel_err = normal.dot(&fluid.velocities[*particle_id]); - if vel_err > na::zero::() { + if vel_err > na::zero::() { fluid.velocities[*particle_id] -= *normal * vel_err; } @@ -237,13 +347,13 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { } let velocity = - body.map(|b| b.velocity_at_point(&proj.point)); + body.map(|b| b.velocity_at_point(proj.point)); - boundary - .velocities - .push(velocity.unwrap_or(Vector::zeros())); - boundary.positions.push(proj.point); - boundary.volumes.push(na::zero::()); + boundary.velocities.push(rapier_to_salva_vector( + velocity.unwrap_or(RapierVector::ZERO), + )); + boundary.positions.push(rapier_to_salva_vector(proj.point)); + boundary.volumes.push(na::zero::()); coupling.features.push(feature); } } @@ -260,7 +370,12 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { } } - fn transmit_forces(&mut self, timestep: &TimestepManager, boundaries: &BoundarySet) { + fn transmit_forces( + &mut self, + timestep: &TimestepManager, + boundaries: &BoundarySet, + force_coefficient: math::Real, + ) { for (collider, coupling) in &self.coupling.entries { if let (Some(collider), Some(boundary)) = ( self.colliders.get(*collider), @@ -277,7 +392,13 @@ impl<'a> CouplingManager for ColliderCouplingManager<'a> { for (pos, force) in boundary.positions.iter().zip(forces.iter().cloned()) { - body.apply_impulse_at_point(force * timestep.dt(), *pos, true) + body.apply_impulse_at_point( + salva_to_rapier_vector( + force * timestep.dt() * force_coefficient, + ), + salva_to_rapier_vector(*pos), + true, + ) } } } diff --git a/src/integrations/rapier/harness_plugin.rs b/src/integrations/rapier/harness_plugin.rs index 8f6c5a4..353e640 100644 --- a/src/integrations/rapier/harness_plugin.rs +++ b/src/integrations/rapier/harness_plugin.rs @@ -57,7 +57,7 @@ impl HarnessPlugin for FluidsHarnessPlugin { } fn step(&mut self, physics: &mut PhysicsState, _run_state: &RunState) { - let step_time = instant::now(); + //let step_time = instant::now(); let dt = physics.integration_parameters.dt; self.fluids_pipeline.step( &physics.gravity, @@ -66,7 +66,7 @@ impl HarnessPlugin for FluidsHarnessPlugin { &mut physics.bodies, ); - self.step_time = instant::now() - step_time; + //self.step_time = instant::now() - step_time; } fn profiling_string(&self) -> String { diff --git a/src/integrations/rapier/mod.rs b/src/integrations/rapier/mod.rs index 7ffd0a5..4e354be 100644 --- a/src/integrations/rapier/mod.rs +++ b/src/integrations/rapier/mod.rs @@ -1,5 +1,6 @@ //! Two-way coupling with the Rapier physics engine. +pub use crate::solver::DfsphParameters; pub use fluids_pipeline::{ ColliderCouplingManager, ColliderCouplingSet, ColliderSampling, FluidsPipeline, }; diff --git a/src/integrations/rapier/testbed_plugin.rs b/src/integrations/rapier/testbed_plugin.rs index 8699dd1..76a2679 100644 --- a/src/integrations/rapier/testbed_plugin.rs +++ b/src/integrations/rapier/testbed_plugin.rs @@ -1,21 +1,21 @@ -use crate::math::{Isometry, Point, Real, Rotation, Translation, Vector}; +use crate::math::{Real, Vector}; use crate::object::{BoundaryHandle, FluidHandle}; -use bevy::math::Quat; -use bevy::prelude::{Assets, Commands, Mesh, Query, Transform}; -use bevy_egui::{egui::ComboBox, egui::Window, EguiContexts}; +#[cfg(feature = "dim2")] +use kiss3d::prelude::Vec2; #[cfg(feature = "dim3")] -use na::Quaternion; -use na::{Point3, Vector3}; -use parry::shape::SharedShape; +use kiss3d::prelude::Vec3; +use kiss3d::{color::Color, window::Window}; +use na::Vector3; use rapier_testbed::{ - harness::Harness, objects::node::EntityWithGraphics, BevyMaterial, GraphicsManager, - PhysicsState, TestbedPlugin, + egui, harness::Harness, settings::ExampleSettings, GraphicsManager, PhysicsState, Testbed, + TestbedPlugin, }; -use crate::integrations::rapier::FluidsPipeline; +use crate::integrations::rapier::{DfsphParameters, FluidsPipeline}; +use std::cell::RefCell; use std::collections::HashMap; +use std::rc::Rc; -//FIXME: handle this with macros, or use bevy-inspectable-egui pub const FLUIDS_RENDERING_MAP: [(&str, FluidsRenderingMode); 3] = [ ("Static", FluidsRenderingMode::StaticColor), ( @@ -25,13 +25,6 @@ pub const FLUIDS_RENDERING_MAP: [(&str, FluidsRenderingMode); 3] = [ max: 50.0, }, ), - // ( - // "Velocity Color & Opacity", - // FluidsRenderingMode::VelocityColorOpacity { - // min: 0.0, - // max: 50.0, - // }, - // ), ( "Velocity Arrows", FluidsRenderingMode::VelocityArrows { @@ -39,9 +32,16 @@ pub const FLUIDS_RENDERING_MAP: [(&str, FluidsRenderingMode); 3] = [ max: 50.0, }, ), - // ("Acceleration Arrows", FluidsRenderingMode::AccelerationArrows), ]; +const SETTINGS_RENDER_BOUNDARIES: &str = "Fluid render boundaries"; +const SETTINGS_RENDERING_MODE: &str = "Fluid rendering mode"; +const SETTINGS_MAX_PRESSURE_ITER: &str = "Fluid max pressure iter"; +const SETTINGS_MAX_DIVERGENCE_ITER: &str = "Fluid max divergence iter"; +const SETTINGS_MAX_DENSITY_ERROR: &str = "Fluid max density error"; +const SETTINGS_MAX_DIVERGENCE_ERROR: &str = "Fluid max divergence error"; +const SETTINGS_BOUNDARY_FORCE_COEFFICIENT: &str = "Fluid boundary force"; + /// How the fluids should be rendered by the testbed. #[derive(Copy, Clone, Debug, PartialEq)] pub enum FluidsRenderingMode { @@ -50,46 +50,39 @@ pub enum FluidsRenderingMode { /// Use a red taint the closer to `max` the velocity is. VelocityColor { /// Fluids with a velocity smaller than this will not have any red taint. - min: f32, + min: Real, /// Fluids with a velocity greater than this will be completely red. - max: f32, + max: Real, }, - // /// Use a red taint the closer to `max` the velocity is, with opacity, low velocity is more transparent - // VelocityColorOpacity { - // /// Fluids with a velocity smaller than this will not have any red taint. - // min: f32, - // /// Fluids with a velocity greater than this will be completely red. - // max: f32, - // }, - /// Show particles as arrows indicating the velocity + /// Show particles as arrows indicating the velocity. VelocityArrows { /// Fluids with a velocity smaller than this will not have any red taint. - min: f32, + min: Real, /// Fluids with a velocity greater than this will be completely red. - max: f32, + max: Real, }, } /// A user-defined callback executed at each frame. pub type FluidCallback = Box; -/// A plugin for rendering fluids with the Rapier testbed. +/// A plugin for stepping fluids inside the Rapier testbed. pub struct FluidsTestbedPlugin { - /// Whether to render the boundary particles + /// Whether to render the boundary particles. pub render_boundary_particles: bool, - /// Rendering mode of fluid particles + /// Rendering mode of fluid particles. pub fluids_rendering_mode: FluidsRenderingMode, callbacks: Vec, step_time: f64, fluids_pipeline: FluidsPipeline, - f2sn: HashMap>, - boundary2sn: HashMap>, - f2color: HashMap>, - ground_color: Point3, - default_fluid_color: Point3, - queue_graphics_reset: bool, + dfsph_parameters: DfsphParameters, + f2color: HashMap>, + boundary2color: HashMap>, + default_fluid_color: Vector3, } +struct SharedFluidsTestbedPlugin(Rc>); + impl FluidsTestbedPlugin { /// Initializes the plugin. pub fn new() -> Self { @@ -99,12 +92,10 @@ impl FluidsTestbedPlugin { step_time: 0.0, callbacks: Vec::new(), fluids_pipeline: FluidsPipeline::new(0.025, 2.0), - f2sn: HashMap::new(), - boundary2sn: HashMap::new(), + dfsph_parameters: DfsphParameters::default(), f2color: HashMap::new(), - ground_color: Point3::new(0.5, 0.5, 0.5), - default_fluid_color: Point3::new(0.0, 0.0, 0.5), - queue_graphics_reset: false, + boundary2color: HashMap::new(), + default_fluid_color: Vector3::new(0.0, 0.0, 0.5), } } @@ -113,14 +104,35 @@ impl FluidsTestbedPlugin { self.callbacks.push(Box::new(f)) } + /// Adds this plugin to the Rapier testbed and registers fluid rendering. + pub fn add_to_testbed(self, testbed: &mut Testbed) { + let plugin = Rc::new(RefCell::new(self)); + let renderer = Rc::clone(&plugin); + + testbed.add_callback(move |graphics, _physics, _events, _run_state| { + if let Some(graphics) = graphics { + let mut renderer = renderer.borrow_mut(); + + if let Some(settings) = graphics.settings.as_deref_mut() { + renderer.update_from_settings(settings); + } + + renderer.draw_fluids(graphics.window); + } + }); + + testbed.add_plugin(SharedFluidsTestbedPlugin(plugin)); + } + /// Sets the fluids pipeline used by the testbed. pub fn set_pipeline(&mut self, fluids_pipeline: FluidsPipeline) { self.fluids_pipeline = fluids_pipeline; self.fluids_pipeline.liquid_world.counters.enable(); + self.refresh_dfsph_parameters(); } /// Sets the color used to render the specified fluid. - pub fn set_fluid_color(&mut self, fluid: FluidHandle, color: Point3) { + pub fn set_fluid_color(&mut self, fluid: FluidHandle, color: Vector3) { let _ = self.f2color.insert(fluid, color); } @@ -134,203 +146,243 @@ impl FluidsTestbedPlugin { self.render_boundary_particles = enabled; } - // TODO: pass velocity & acceleration vectors in - fn add_particle_graphics( - &self, - particle: &Point, - particle_radius: f32, - graphics: &mut GraphicsManager, - commands: &mut Commands, - meshes: &mut Assets, - materials: &mut Assets, - _components: &mut Query<&mut Transform>, - _harness: &mut Harness, - color: &Point3, - force_shape: Option, - ) -> Vec { - let shape = if let Some(shape) = force_shape { - shape - } else { - match self.fluids_rendering_mode { - #[cfg(feature = "dim3")] - FluidsRenderingMode::VelocityArrows { .. } => { - SharedShape::cone(particle_radius, particle_radius / 4.) - } - // #[cfg(feature = "dim2")] - //FIXME: This doesn't work, it is caused by either not being in prefab_meshes, or the shape_type not being supported.. somewhere - // FluidsRenderingMode::VelocityArrows { .. } => SharedShape::triangle( - // Point::new(0., particle_radius), - // Point::new(particle_radius * 0.4, -particle_radius * 0.8), - // Point::new(-particle_radius * 0.4, -particle_radius * 0.8), - // ), - _ => SharedShape::ball(particle_radius), + fn color(color: Vector3) -> Color { + Color::new(color.x as f32, color.y as f32, color.z as f32, 1.0) + } + + fn fluid_color(&self, handle: FluidHandle, velocity: &Vector) -> Color { + let base = *self + .f2color + .get(&handle) + .unwrap_or(&self.default_fluid_color); + + match self.fluids_rendering_mode { + FluidsRenderingMode::StaticColor => Self::color(base), + FluidsRenderingMode::VelocityColor { min, max } + | FluidsRenderingMode::VelocityArrows { min, max } => { + Self::color(Self::velocity_tinted_color(base, velocity, min, max)) } - }; + } + } + + fn velocity_tinted_color( + base: Vector3, + velocity: &Vector, + min: Real, + max: Real, + ) -> Vector3 { + let range = max - min; + + if range <= na::zero() { + return base; + } + + let factor = ((velocity.norm() - min) / range) + .max(na::zero()) + .min(na::one()); + let inv_factor = na::one::() - factor; + + Vector3::new( + base.x * inv_factor + factor, + base.y * inv_factor, + base.z * inv_factor, + ) + } + + fn particle_size(radius: Real) -> f32 { + (radius as f32 * 300.0).clamp(2.5, 8.0) + } + + fn refresh_dfsph_parameters(&mut self) { + if let Some(parameters) = self.fluids_pipeline.dfsph_parameters() { + self.dfsph_parameters = parameters; + } + } + + fn rendering_mode_index(&self) -> usize { + FLUIDS_RENDERING_MAP + .iter() + .position(|(_, mode)| *mode == self.fluids_rendering_mode) + .unwrap_or(0) + } + + fn get_clamped_u32_setting( + settings: &mut ExampleSettings, + key: &'static str, + default: u32, + min: u32, + max: u32, + ) -> u32 { + let value = settings + .get_or_set_u32(key, default.clamp(min, max), min..=max) + .clamp(min, max); + settings.set_u32(key, value, min..=max); + value + } - let mut shapes = Vec::new(); - let isometry = - Isometry::from_parts(Translation::from(particle.coords), Rotation::identity()); - graphics.add_shape( - commands, - meshes, - materials, - None, - &*shape, - false, - &isometry, - &Isometry::identity(), - *color, - &mut shapes, + fn update_from_settings(&mut self, settings: &mut ExampleSettings) { + self.render_boundary_particles = + settings.get_or_set_bool(SETTINGS_RENDER_BOUNDARIES, self.render_boundary_particles); + + let rendering_options = FLUIDS_RENDERING_MAP + .iter() + .map(|(name, _)| (*name).to_string()) + .collect(); + let rendering_mode = settings.get_or_set_string( + SETTINGS_RENDERING_MODE, + self.rendering_mode_index(), + rendering_options, ); - shapes + + if let Some((_, mode)) = FLUIDS_RENDERING_MAP.get(rendering_mode) { + self.fluids_rendering_mode = *mode; + } + + let max_pressure_iter = Self::get_clamped_u32_setting( + settings, + SETTINGS_MAX_PRESSURE_ITER, + self.dfsph_parameters.max_pressure_iter as u32, + 1, + 80, + ) as usize; + let max_divergence_iter = Self::get_clamped_u32_setting( + settings, + SETTINGS_MAX_DIVERGENCE_ITER, + self.dfsph_parameters.max_divergence_iter as u32, + 1, + 80, + ) as usize; + let max_density_error = settings.get_or_set_f32( + SETTINGS_MAX_DENSITY_ERROR, + self.dfsph_parameters.max_density_error as f32, + 0.0..=0.5, + ); + let max_divergence_error = settings.get_or_set_f32( + SETTINGS_MAX_DIVERGENCE_ERROR, + self.dfsph_parameters.max_divergence_error as f32, + 0.0..=2.0, + ); + + let boundary_force = settings.get_or_set_f32( + SETTINGS_BOUNDARY_FORCE_COEFFICIENT, + self.fluids_pipeline.liquid_world.boundary_force_coefficient as f32, + 0.0..=1.0, + ); + self.fluids_pipeline.liquid_world.boundary_force_coefficient = + na::convert::<_, Real>(boundary_force); + + let dfsph_parameters = DfsphParameters { + min_pressure_iter: self.dfsph_parameters.min_pressure_iter, + max_pressure_iter, + max_density_error: na::convert::<_, Real>(max_density_error), + min_divergence_iter: self.dfsph_parameters.min_divergence_iter, + max_divergence_iter, + max_divergence_error: na::convert::<_, Real>(max_divergence_error), + }; + + if dfsph_parameters != self.dfsph_parameters { + self.dfsph_parameters = dfsph_parameters; + let _ = self.fluids_pipeline.set_dfsph_parameters(dfsph_parameters); + } } - fn lerp_velocity( - velocity: Vector, - start: Vector3, - min: f32, - max: f32, - ) -> Vector3 { - let end = Vector3::new(1.0, 0.0, 0.0); - let vel: Vector = na::convert_unchecked(velocity); - let vel: Vector = na::convert(vel); - let t = (vel.norm() - min) / (max - min); - start.lerp(&end, na::clamp(t, 0.0, 1.0)) + #[cfg(feature = "dim2")] + fn point(point: &Vector) -> Vec2 { + Vec2::new(point.x as f32, point.y as f32) } -} -impl TestbedPlugin for FluidsTestbedPlugin { - fn init_plugin(&mut self) { - // TODO: decide if anything needs to be changed + #[cfg(feature = "dim3")] + fn point(point: &Vector) -> Vec3 { + Vec3::new(point.x as f32, point.y as f32, point.z as f32) } - fn init_graphics( - &mut self, - graphics: &mut GraphicsManager, - commands: &mut Commands, - meshes: &mut Assets, - materials: &mut Assets, - components: &mut Query<&mut Transform>, - harness: &mut Harness, + #[cfg(feature = "dim2")] + fn draw_particle(window: &mut Window, point: &Vector, color: Color, size: f32) { + window.draw_point_2d(Self::point(point), color, size); + } + + #[cfg(feature = "dim3")] + fn draw_particle(window: &mut Window, point: &Vector, color: Color, size: f32) { + window.draw_point(Self::point(point), color, size); + } + + #[cfg(feature = "dim2")] + fn draw_velocity( + window: &mut Window, + point: &Vector, + velocity: &Vector, + color: Color, ) { + let end = *point + *velocity * na::convert::<_, Real>(0.02); + window.draw_line_2d(Self::point(point), Self::point(&end), color, 1.5); + } + + #[cfg(feature = "dim3")] + fn draw_velocity( + window: &mut Window, + point: &Vector, + velocity: &Vector, + color: Color, + ) { + let end = *point + *velocity * na::convert::<_, Real>(0.02); + window.draw_line(Self::point(point), Self::point(&end), color, 1.5, false); + } + + fn draw_fluids(&self, window: &mut Window) { + let draw_velocities = matches!( + self.fluids_rendering_mode, + FluidsRenderingMode::VelocityArrows { .. } + ); + for (handle, fluid) in self.fluids_pipeline.liquid_world.fluids().iter() { - let _ = self - .f2sn - .insert(handle, Vec::with_capacity(fluid.positions.len())); - - let color = *self - .f2color - .entry(handle) - .or_insert_with(|| self.default_fluid_color); - - for particle in &fluid.positions { - let ent = self.add_particle_graphics( - particle, - fluid.particle_radius(), - graphics, - commands, - meshes, - materials, - components, - harness, - &color, - None, - ); - if let Some(entities) = self.f2sn.get_mut(&handle) { - entities.extend(ent); + let size = Self::particle_size(fluid.particle_radius()); + + for (point, velocity) in fluid.positions.iter().zip(fluid.velocities.iter()) { + let color = self.fluid_color(handle, velocity); + Self::draw_particle(window, point, color, size); + + if draw_velocities && velocity.norm_squared() > na::zero() { + Self::draw_velocity(window, point, velocity, color); } } } - let particle_radius = self.fluids_pipeline.liquid_world.particle_radius(); - - // FIXME: There is currently no way to get the collider pose from this function if self.render_boundary_particles { + let default_color = Vector3::repeat(na::convert::<_, Real>(0.5)); + let size = Self::particle_size(self.fluids_pipeline.liquid_world.particle_radius()); + for (handle, boundary) in self.fluids_pipeline.liquid_world.boundaries().iter() { - let _ = self - .boundary2sn - .insert(handle, Vec::with_capacity(boundary.num_particles())); - let color = self.ground_color; - - for (_, cce) in &self.fluids_pipeline.coupling.entries { - if cce.boundary == handle { - match &cce.sampling_method { - crate::integrations::rapier::ColliderSampling::StaticSampling(particles) => { - for particle in particles { - let ent = self.add_particle_graphics( - particle, - particle_radius, - graphics, - commands, - meshes, - materials, - components, - harness, - &color, - Some(SharedShape::ball(particle_radius)) - ); - if let Some(entities) = self.boundary2sn.get_mut(&handle) { - entities.extend(ent); - } - } - }, - crate::integrations::rapier::ColliderSampling::DynamicContactSampling => { - // TODO: ??? - }, - } - } - } - } - } - } + let color = + Self::color(*self.boundary2color.get(&handle).unwrap_or(&default_color)); - fn clear_graphics(&mut self, _graphics: &mut GraphicsManager, commands: &mut Commands) { - for (handle, _) in self.fluids_pipeline.liquid_world.fluids().iter() { - if let Some(entities) = self.f2sn.get_mut(&handle) { - for entity in entities { - entity.despawn(commands); + for point in &boundary.positions { + Self::draw_particle(window, point, color, size); } } } + } +} - for (handle, _) in self.fluids_pipeline.liquid_world.boundaries().iter() { - if let Some(entities) = self.boundary2sn.get_mut(&handle) { - for entity in entities { - entity.despawn(commands); - } - } - } +impl TestbedPlugin for FluidsTestbedPlugin { + fn init_plugin(&mut self) {} - self.f2sn.clear(); - self.boundary2sn.clear(); + fn init_graphics( + &mut self, + _graphics: &mut GraphicsManager, + _window: &mut Window, + _harness: &mut Harness, + ) { } - fn run_callbacks(&mut self, harness: &mut Harness) { - // FIXME: salva should be able to keep a list of indices that were added & removed in this step - // at the moment we just clear & initialize the grahics when fluids_lengths changes - let fluid_lengths: Vec<(FluidHandle, usize)> = self - .fluids_pipeline - .liquid_world - .fluids() - .iter() - .map(|(h, f)| (h, f.positions.len())) - .collect(); + fn clear_graphics(&mut self, _graphics: &mut GraphicsManager, _window: &mut Window) {} + fn run_callbacks(&mut self, harness: &mut Harness) { for f in &mut self.callbacks { f(harness, &mut self.fluids_pipeline) } - - for (h, fl) in fluid_lengths { - if let Some(fluid) = self.fluids_pipeline.liquid_world.fluids().get(h) { - if fluid.positions.len() != fl || fluid.num_deleted_particles() > 0 { - self.queue_graphics_reset = true; - } - } - } } fn step(&mut self, physics: &mut PhysicsState) { - let step_time = instant::now(); let dt = physics.integration_parameters.dt; self.fluids_pipeline.step( &physics.gravity, @@ -338,174 +390,93 @@ impl TestbedPlugin for FluidsTestbedPlugin { &physics.colliders, &mut physics.bodies, ); - - self.step_time = instant::now() - step_time; } fn draw( + &mut self, + _graphics: &mut GraphicsManager, + window: &mut Window, + _harness: &mut Harness, + ) { + self.draw_fluids(window); + } + + fn update_ui( + &mut self, + ui_context: &egui::Context, + _harness: &mut Harness, + _graphics: &mut GraphicsManager, + _window: &mut Window, + ) { + let _ = egui::Window::new("Fluids").show(ui_context, |ui| { + let _ = ui.checkbox( + &mut self.render_boundary_particles, + "Render boundary particles", + ); + + let selected = FLUIDS_RENDERING_MAP + .iter() + .find_map(|(name, mode)| (*mode == self.fluids_rendering_mode).then_some(*name)) + .unwrap_or("Custom"); + + let _ = egui::ComboBox::from_label("Rendering mode") + .selected_text(selected) + .show_ui(ui, |ui| { + for (name, mode) in FLUIDS_RENDERING_MAP { + let _ = ui.selectable_value(&mut self.fluids_rendering_mode, mode, name); + } + }); + }); + } + + fn profiling_string(&self) -> String { + format!("Fluids: {:.2}ms", self.step_time) + } +} + +impl TestbedPlugin for SharedFluidsTestbedPlugin { + fn init_plugin(&mut self) { + self.0.borrow_mut().init_plugin(); + } + + fn init_graphics( &mut self, graphics: &mut GraphicsManager, - commands: &mut Commands, - meshes: &mut Assets, - materials: &mut Assets, - components: &mut Query<&mut Transform>, + window: &mut Window, harness: &mut Harness, ) { - if self.queue_graphics_reset { - self.clear_graphics(graphics, commands); - self.init_graphics(graphics, commands, meshes, materials, components, harness); - self.queue_graphics_reset = false; - } + self.0.borrow_mut().init_graphics(graphics, window, harness); + } - let (mut min, mut max) = (f32::MAX, f32::MIN); - for (handle, fluid) in self.fluids_pipeline.liquid_world.fluids().iter() { - if let Some(entities) = self.f2sn.get_mut(&handle) { - for (idx, particle) in fluid.positions.iter().enumerate() { - let velocity = Vector::from(fluid.velocities[idx]); - let magnitude = velocity.magnitude(); - min = min.min(magnitude); - max = max.max(magnitude); - if let Some(entity) = entities.get_mut(idx) { - if let Ok(mut pos) = components.get_mut(entity.entity) { - { - pos.translation.x = particle.x; - pos.translation.y = particle.y; - #[cfg(feature = "dim3")] - { - pos.translation.z = particle.z; - - if let FluidsRenderingMode::VelocityArrows { .. } = - self.fluids_rendering_mode - { - let cone_paxis: Quaternion = - Quaternion::from_vector(-Vector3::y().to_homogeneous()); - let vr = Quaternion::from_vector( - velocity.normalize().to_homogeneous(), - ); - let rotation = (vr - cone_paxis).normalize(); - - pos.rotation = Quat::from_xyzw( - rotation.i, rotation.j, rotation.k, rotation.w, - ); - } - } - #[cfg(feature = "dim2")] - { - let norm = velocity.normalize(); - let hyp = (norm.x * norm.x + norm.y * norm.y).sqrt(); - let angle = 2. * (norm.y / (norm.x + hyp)).atan(); - pos.rotation = Quat::from_rotation_z(angle); - } - } - } - - if let Some(color) = self.f2color.get(&handle) { - match self.fluids_rendering_mode { - FluidsRenderingMode::VelocityColor { min, max } => { - let lerp = Self::lerp_velocity( - fluid.velocities[idx], - color.coords, - min, - max, - ); - entity - .set_color(materials, Point3::new(lerp.x, lerp.y, lerp.z)); - } - FluidsRenderingMode::VelocityArrows { min, max } => { - let lerp = Self::lerp_velocity( - fluid.velocities[idx], - color.coords, - min, - max, - ); - entity - .set_color(materials, Point3::new(lerp.x, lerp.y, lerp.z)); - } - // FIXME: rapier needs to be updated to respect opacity - // FluidsRenderingMode::VelocityColorOpacity { min, max } => { - // let lerp = lerp_velocity( - // fluid.velocities[idx], - // color.coords, - // min, - // max, - // ); - // entity.opacity = lerp.magnitude(); - // entity.set_color( - // _materials, - // Point3::new(lerp.x, lerp.y, lerp.z), - // ); - // } - // FluidsRenderingMode::VelocityArrows => {} - _ => { - entity.opacity = 1.0; - entity.set_color(materials, *color); - } - } - } - entity.update(&harness.physics.colliders, components, &graphics.gfx_shift) - } - } - } - } + fn clear_graphics(&mut self, graphics: &mut GraphicsManager, window: &mut Window) { + self.0.borrow_mut().clear_graphics(graphics, window); + } + + fn run_callbacks(&mut self, harness: &mut Harness) { + self.0.borrow_mut().run_callbacks(harness); + } + + fn step(&mut self, physics: &mut PhysicsState) { + self.0.borrow_mut().step(physics); + } + + fn draw(&mut self, graphics: &mut GraphicsManager, window: &mut Window, harness: &mut Harness) { + self.0.borrow_mut().draw(graphics, window, harness); } fn update_ui( &mut self, - ui_context: &EguiContexts, + ui_context: &egui::Context, harness: &mut Harness, graphics: &mut GraphicsManager, - commands: &mut Commands, - meshes: &mut Assets, - materials: &mut Assets, - components: &mut Query<&mut Transform>, + window: &mut Window, ) { - fn get_rendering_mode_index(rendering_mode: FluidsRenderingMode) -> usize { - FLUIDS_RENDERING_MAP - .iter() - .enumerate() - .find(|(_, mode)| rendering_mode == mode.1) - .map(|(idx, _)| idx) - .unwrap_or(0) - } - - let _ = Window::new("Fluid Parameters") - .min_height(200.0) - .show(ui_context.ctx(), |ui| { - let mut changed = false; - - let _ = ComboBox::from_label("Rendering Mode") - .width(150.0) - .selected_text( - FLUIDS_RENDERING_MAP[get_rendering_mode_index(self.fluids_rendering_mode)] - .0, - ) - .show_ui(ui, |ui| { - for (_, (name, mode)) in FLUIDS_RENDERING_MAP.iter().enumerate() { - changed = ui - .selectable_value(&mut self.fluids_rendering_mode, *mode, *name) - .changed() - || changed; - } - }); - - if changed { - // FIXME: not too sure what to do here for color - // let fluid_handle = self - // .fluids_pipeline - // .liquid_world - // .fluids() - // .iter() - // .next() - // .unwrap() - // .0; - self.clear_graphics(graphics, commands); - // let color = self.f2color[&fluid_handle].clone(); - self.init_graphics(graphics, commands, meshes, materials, components, harness) - } - }); + self.0 + .borrow_mut() + .update_ui(ui_context, harness, graphics, window); } fn profiling_string(&self) -> String { - format!("Fluids: {:.2}ms", self.step_time) + self.0.borrow().profiling_string() } } diff --git a/src/kernel/kernel.rs b/src/kernel/kernel.rs index 34d430f..486597b 100644 --- a/src/kernel/kernel.rs +++ b/src/kernel/kernel.rs @@ -1,4 +1,4 @@ -use crate::math::{Point, Real, Vector}; +use crate::math::{Real, Vector}; use approx::AbsDiffEq; use na::Unit; @@ -24,17 +24,17 @@ pub trait Kernel: Send + Sync { } /// Evaluate the kernel for the vector equal to `p1 - p2`. - fn points_apply(p1: &Point, p2: &Point, h: Real) -> Real { + fn points_apply(p1: &Vector, p2: &Vector, h: Real) -> Real { Self::apply(p1 - p2, h) } /// Differential wrt. the coordinates of `p1`. - fn points_apply_diff1(p1: &Point, p2: &Point, h: Real) -> Vector { + fn points_apply_diff1(p1: &Vector, p2: &Vector, h: Real) -> Vector { Self::apply_diff(p1 - p2, h) } /// Differential wrt. the coordinates of `p2`. - fn points_apply_diff2(p1: &Point, p2: &Point, h: Real) -> Vector { + fn points_apply_diff2(p1: &Vector, p2: &Vector, h: Real) -> Vector { -Self::apply_diff(p1 - p2, h) } } diff --git a/src/lib.rs b/src/lib.rs index 650ee92..7a44d56 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -26,7 +26,7 @@ The name of this library is inspired from the famous surrealist artist `Salvador #![deny(non_camel_case_types)] #![deny(unused_parens)] #![deny(non_upper_case_globals)] -#![deny(unused_qualifications)] +#![allow(unused_qualifications)] #![warn(missing_docs)] // FIXME: deny this #![deny(unused_results)] #![allow(type_alias_bounds)] @@ -37,18 +37,30 @@ The name of this library is inspired from the famous surrealist artist `Salvador extern crate nalgebra as na; extern crate num_traits as num; -#[cfg(all(feature = "dim2", feature = "parry"))] +#[cfg(all(feature = "dim2", feature = "f32", feature = "parry"))] pub extern crate parry2d as parry; -#[cfg(all(feature = "dim3", feature = "parry"))] +#[cfg(all(feature = "dim2", feature = "f64", feature = "parry"))] +pub extern crate parry2d_f64 as parry; +#[cfg(all(feature = "dim3", feature = "f32", feature = "parry"))] pub extern crate parry3d as parry; -#[cfg(all(feature = "dim2", feature = "rapier"))] +#[cfg(all(feature = "dim3", feature = "f64", feature = "parry"))] +pub extern crate parry3d_f64 as parry; +#[cfg(all(feature = "dim2", feature = "f32", feature = "rapier"))] pub extern crate rapier2d as rapier; -#[cfg(all(feature = "dim3", feature = "rapier"))] +#[cfg(all(feature = "dim2", feature = "f64", feature = "rapier"))] +pub extern crate rapier2d_f64 as rapier; +#[cfg(all(feature = "dim3", feature = "f32", feature = "rapier"))] pub extern crate rapier3d as rapier; -#[cfg(all(feature = "dim2", feature = "rapier-testbed"))] -extern crate rapier_testbed2d as rapier_testbed; -#[cfg(all(feature = "dim3", feature = "rapier-testbed"))] -extern crate rapier_testbed3d as rapier_testbed; +#[cfg(all(feature = "dim3", feature = "f64", feature = "rapier"))] +pub extern crate rapier3d_f64 as rapier; +#[cfg(all(feature = "dim2", feature = "f32", feature = "rapier-testbed"))] +pub extern crate rapier_testbed2d as rapier_testbed; +#[cfg(all(feature = "dim2", feature = "f64", feature = "rapier-testbed"))] +pub extern crate rapier_testbed2d_f64 as rapier_testbed; +#[cfg(all(feature = "dim3", feature = "f32", feature = "rapier-testbed"))] +pub extern crate rapier_testbed3d as rapier_testbed; +#[cfg(all(feature = "dim3", feature = "f64", feature = "rapier-testbed"))] +pub extern crate rapier_testbed3d_f64 as rapier_testbed; macro_rules! par_iter { ($t: expr) => {{ @@ -103,7 +115,7 @@ pub use crate::timestep_manager::TimestepManager; #[cfg(feature = "dim3")] pub mod math { use na::{ - Isometry3, Matrix3, Matrix6, Matrix6xX, MatrixView6xX, MatrixViewMut6xX, Point3, Rotation3, + Isometry3, Matrix3, Matrix6, Matrix6xX, MatrixView6xX, MatrixViewMut6xX, Rotation3, Translation3, UnitQuaternion, Vector3, Vector6, U3, U6, }; @@ -114,8 +126,12 @@ pub mod math { /// The maximum number of possible translations of a rigid body. pub const DIM: usize = 3; + #[cfg(all(feature = "f32"))] /// The scalar type. pub type Real = f32; + #[cfg(all(feature = "f64"))] + /// The scalar type. + pub type Real = f64; /// The dimension of the ambient space. pub type Dim = U3; @@ -126,9 +142,6 @@ pub mod math { /// The dimension of the rotations. pub type AngularDim = U3; - /// The point type. - pub type Point = Point3; - /// The angular vector type. pub type AngularVector = Vector3; @@ -184,7 +197,7 @@ pub mod math { #[cfg(feature = "dim2")] pub mod math { use na::{ - Isometry2, Matrix1, Matrix2, Matrix3, Matrix6xX, MatrixView3xX, MatrixViewMut3xX, Point2, + Isometry2, Matrix1, Matrix2, Matrix3, Matrix6xX, MatrixView3xX, MatrixViewMut3xX, Rotation2, RowVector2, Translation2, UnitComplex, Vector1, Vector2, Vector3, U1, U2, U3, }; @@ -195,8 +208,12 @@ pub mod math { /// The maximum number of possible translations of a rigid body. pub const DIM: usize = 2; + #[cfg(all(feature = "f32"))] /// The scalar type. pub type Real = f32; + #[cfg(all(feature = "f64"))] + /// The scalar type. + pub type Real = f64; /// The dimension of the ambient space. pub type Dim = U2; @@ -207,9 +224,6 @@ pub mod math { /// The dimension of a spatial vector. pub type SpatialDim = U3; - /// The point type. - pub type Point = Point2; - /// The vector type with dimension `SpatialDim × 1`. pub type SpatialVector = Vector3; diff --git a/src/liquid_world.rs b/src/liquid_world.rs index 7ebe81a..610e3af 100644 --- a/src/liquid_world.rs +++ b/src/liquid_world.rs @@ -4,7 +4,7 @@ use crate::geometry::{self, ContactManager, HGrid, HGridEntry}; use crate::math::{Real, Vector}; use crate::object::{Boundary, BoundaryHandle, BoundarySet}; use crate::object::{Fluid, FluidHandle, FluidSet}; -use crate::solver::PressureSolver; +use crate::solver::{DfsphParameters, PressureSolver}; use crate::TimestepManager; #[cfg(feature = "parry")] use { @@ -26,6 +26,10 @@ pub struct LiquidWorld { contact_manager: ContactManager, timestep_manager: TimestepManager, hgrid: HGrid, + /// Coefficient applied when transmitting forces from fluids to coupled rigid bodies. + /// A value of 1.0 applies full force, 0.5 applies half force, etc. + /// Default is 1.0. + pub boundary_force_coefficient: Real, } impl LiquidWorld { @@ -36,10 +40,13 @@ impl LiquidWorld { /// - `particle_radius`: the radius of every particle on this world. /// - `smoothing_factor`: the smoothing factor used to compute the SPH kernel radius. /// The kernel radius will be computed as `particle_radius * smoothing_factor * 2.0. + /// - `boundary_force_coefficient`: coefficient applied when transmitting forces from fluids to boundaries. + /// Default is 1.0 (full force). Use lower values (e.g., 0.5) to reduce fluid influence on rigid bodies. pub fn new( solver: impl PressureSolver + Send + Sync + 'static, particle_radius: Real, smoothing_factor: Real, + boundary_force_coefficient: Real, ) -> Self { let h = particle_radius * smoothing_factor * na::convert::<_, Real>(2.0); Self { @@ -53,6 +60,7 @@ impl LiquidWorld { contact_manager: ContactManager::new(), timestep_manager: TimestepManager::new(particle_radius), hgrid: HGrid::new(h), + boundary_force_coefficient, } } @@ -63,6 +71,16 @@ impl LiquidWorld { self.step_with_coupling(dt, gravity, &mut ()) } + /// Current DFSPH tuning parameters, if this world uses DFSPH. + pub fn dfsph_parameters(&self) -> Option { + self.solver.dfsph_parameters() + } + + /// Set DFSPH tuning parameters if this world uses DFSPH. + pub fn set_dfsph_parameters(&mut self, parameters: DfsphParameters) -> bool { + self.solver.set_dfsph_parameters(parameters) + } + /// Advances the simulation by `dt` seconds, taking into account coupling with an external rigid-body engine. pub fn step_with_coupling( &mut self, @@ -143,7 +161,11 @@ impl LiquidWorld { self.boundaries.as_slice(), ); - coupling.transmit_forces(&self.timestep_manager, &self.boundaries); + coupling.transmit_forces( + &self.timestep_manager, + &self.boundaries, + self.boundary_force_coefficient, + ); self.counters.stages.solver_time.pause(); } @@ -213,17 +235,23 @@ impl LiquidWorld { &'a self, aabb: Aabb, ) -> impl Iterator + 'a { + let mins = aabb.mins.into(); + let maxs = aabb.maxs.into(); self.hgrid - .cells_intersecting_aabb(&aabb.mins, &aabb.maxs) + .cells_intersecting_aabb(&mins, &maxs) .flat_map(|e| e.1) .filter_map(move |entry| match entry { HGridEntry::FluidParticle(fid, pid) => { let (fluid, handle) = self.fluids.get_from_contiguous_index(*fid)?; + + // Defensive bounds check for fluids + if *pid >= fluid.positions.len() { + return None; + } + let pt = fluid.positions[*pid]; - // FIXME: use `distance_to_local_point` once it's supported. - let id = &Isometry::identity(); - if aabb.distance_to_point(id, &pt, true) < self.particle_radius { + if aabb.distance_to_local_point(pt.into(), true) < self.particle_radius { Some(ParticleId::FluidParticle(handle, *pid)) } else { None @@ -231,9 +259,14 @@ impl LiquidWorld { } HGridEntry::BoundaryParticle(bid, pid) => { let (boundary, handle) = self.boundaries.get_from_contiguous_index(*bid)?; - let pt = boundary.positions[*pid]; // FIXME: use `distance_to_local_point` once it's supported. - let id = &Isometry::identity(); - if aabb.distance_to_point(id, &pt, true) < self.particle_radius { + + // Defensive bounds check for fluids + if *pid >= boundary.positions.len() { + return None; + } + + let pt = boundary.positions[*pid]; + if aabb.distance_to_local_point(pt.into(), true) < self.particle_radius { Some(ParticleId::BoundaryParticle(handle, *pid)) } else { None @@ -252,16 +285,25 @@ impl LiquidWorld { where S: Shape, { - let aabb = shape.compute_aabb(pos); + let pos = (*pos).into(); + let aabb = shape.compute_aabb(&pos); + let mins = aabb.mins.into(); + let maxs = aabb.maxs.into(); self.hgrid - .cells_intersecting_aabb(&aabb.mins, &aabb.maxs) + .cells_intersecting_aabb(&mins, &maxs) .flat_map(|e| e.1) .filter_map(move |entry| match entry { HGridEntry::FluidParticle(fid, pid) => { let (fluid, handle) = self.fluids.get_from_contiguous_index(*fid)?; + + // Defensive bounds check for fluids + if *pid >= fluid.positions.len() { + return None; + } + let pt = fluid.positions[*pid]; - if shape.distance_to_point(pos, &pt, true) <= self.particle_radius { + if shape.distance_to_point(&pos, pt.into(), true) <= self.particle_radius { Some(ParticleId::FluidParticle(handle, *pid)) } else { None @@ -269,8 +311,14 @@ impl LiquidWorld { } HGridEntry::BoundaryParticle(bid, pid) => { let (boundary, handle) = self.boundaries.get_from_contiguous_index(*bid)?; - let pt = boundary.positions[*pid]; // FIXME: use `distance_to_local_point` once it's supported. - if shape.distance_to_point(pos, &pt, true) <= self.particle_radius { + + // Defensive bounds check for fluids + if *pid >= boundary.positions.len() { + return None; + } + + let pt = boundary.positions[*pid]; + if shape.distance_to_point(&pos, pt.into(), true) <= self.particle_radius { Some(ParticleId::BoundaryParticle(handle, *pid)) } else { None diff --git a/src/object/boundary.rs b/src/object/boundary.rs index 006bf95..4a41d9f 100644 --- a/src/object/boundary.rs +++ b/src/object/boundary.rs @@ -1,4 +1,4 @@ -use crate::math::{Isometry, Point, Real, Vector}; +use crate::math::{Isometry, Real, Vector}; use crate::object::{ContiguousArena, ContiguousArenaIndex}; use std::sync::RwLock; @@ -10,7 +10,7 @@ use super::interaction_groups::InteractionGroups; /// A boundary object is composed of static particles, or of particles coupled with non-fluid bodies. pub struct Boundary { /// The world-space position of the boundary particles. - pub positions: Vec>, + pub positions: Vec>, /// The artificial velocities of each boundary particle. pub velocities: Vec>, /// The volume computed for each boundary particle. @@ -26,7 +26,7 @@ pub struct Boundary { impl Boundary { /// Initialize a boundary object with the given particles. pub fn new( - particle_positions: Vec>, + particle_positions: Vec>, interaction_groups: InteractionGroups, ) -> Self { let num_particles = particle_positions.len(); @@ -53,7 +53,9 @@ impl Boundary { /// Transforms all the particle positions of this boundary by the given isometry. pub fn transform_by(&mut self, pose: &Isometry) { - self.positions.iter_mut().for_each(|p| *p = pose * *p); + self.positions + .iter_mut() + .for_each(|p| *p = pose.translation.vector + pose.rotation * *p); } /// Apply a force `f` to the `i`-th particle of this boundary object. diff --git a/src/object/contiguous_arena.rs b/src/object/contiguous_arena.rs index e7f664f..2fa8f96 100644 --- a/src/object/contiguous_arena.rs +++ b/src/object/contiguous_arena.rs @@ -66,13 +66,13 @@ impl ContiguousArena { #[inline] /// Gets references to all the objects on this set. - pub fn values(&self) -> std::slice::Iter { + pub fn values(&self) -> std::slice::Iter<'_, T> { self.objects.iter() } #[inline] /// Gets mutable references to all the objects on this set. - pub fn values_mut(&mut self) -> std::slice::IterMut { + pub fn values_mut(&mut self) -> std::slice::IterMut<'_, T> { self.objects.iter_mut() } diff --git a/src/object/fluid.rs b/src/object/fluid.rs index 69bb356..dea00fb 100644 --- a/src/object/fluid.rs +++ b/src/object/fluid.rs @@ -1,4 +1,4 @@ -use crate::math::{Isometry, Point, Real, Vector}; +use crate::math::{Isometry, Real, Vector}; use crate::object::{ContiguousArena, ContiguousArenaIndex}; use crate::solver::NonPressureForce; @@ -13,7 +13,7 @@ pub struct Fluid { /// Nonpressure forces this fluid is subject to. pub nonpressure_forces: Vec>, /// The world-space position of the fluid particles. - pub positions: Vec>, + pub positions: Vec>, /// The velocities of the fluid particles. pub velocities: Vec>, /// The accelerations of the fluid particles. @@ -23,7 +23,7 @@ pub struct Fluid { /// The rest density of this fluid. pub density0: Real, /// Mask indicating what particles have been deleted. - deleted_particles: Vec, + pub deleted_particles: Vec, /// Indicates if a bit of the `deleted_particles` mask has been set. num_deleted_particles: usize, /// The particles radius. @@ -38,7 +38,7 @@ impl Fluid { /// /// The particle radius should be the same as the radius used to initialize the liquid world. pub fn new( - particle_positions: Vec>, + particle_positions: Vec>, particle_radius: Real, // XXX: remove this parameter since it is already defined by the liquid world. density0: Real, interaction_groups: InteractionGroups, @@ -125,7 +125,7 @@ impl Fluid { /// If it is not `None`, then it must be a slice with the same length than `positions`. pub fn add_particles( &mut self, - positions: &[Point], + positions: &[Vector], velocities: Option<&[Vector]>, ) { let nparticles = self.positions.len() + positions.len(); @@ -164,7 +164,9 @@ impl Fluid { /// Apply the given transformation to each particle of this fluid. pub fn transform_by(&mut self, t: &Isometry) { - self.positions.iter_mut().for_each(|p| *p = t * *p) + self.positions + .iter_mut() + .for_each(|p| *p = t.translation.vector + t.rotation * *p) } /// The number of particles on this fluid. @@ -176,7 +178,8 @@ impl Fluid { #[cfg(feature = "parry")] pub fn compute_aabb(&self, particle_radius: Real) -> parry::bounding_volume::Aabb { use parry::bounding_volume::{details::local_point_cloud_aabb, BoundingVolume}; - local_point_cloud_aabb(&self.positions).loosened(particle_radius) + local_point_cloud_aabb(self.positions.iter().copied().map(Into::into)) + .loosened(particle_radius) } /// The mass of the `i`-th particle of this fluid. diff --git a/src/object/interaction_groups.rs b/src/object/interaction_groups.rs index 0873168..ab3e68e 100644 --- a/src/object/interaction_groups.rs +++ b/src/object/interaction_groups.rs @@ -7,14 +7,14 @@ /// - The interaction groups memberships. /// - The interaction groups filter. /// -/// An interaction is allowed between two filters `a` and `b` when two conditions -/// are met simultaneously: +/// An interaction is allowed between two filters `a` and `b` when at least one of +/// these conditions is met: /// - The groups membership of `a` has at least one bit set to `1` in common with the groups filter of `b`. /// - The groups membership of `b` has at least one bit set to `1` in common with the groups filter of `a`. /// /// In other words, interactions are allowed between two filter iff. the following condition is met: /// ```ignore -/// (self.memberships & rhs.filter) != 0 && (rhs.memberships & self.filter) != 0 +/// (self.memberships & rhs.filter) != 0 || (rhs.memberships & self.filter) != 0 /// ``` #[derive(Copy, Clone, Debug, Hash, PartialEq, Eq)] #[repr(C)] @@ -58,14 +58,15 @@ impl InteractionGroups { /// Check if interactions should be allowed based on the interaction memberships and filter. /// - /// An interaction is allowed iff. the memberships of `self` contain at least one bit set to 1 in common - /// with the filter of `rhs`, and vice-versa. + /// An interaction is allowed if the memberships of `self` contain at least one bit set to 1 in common + /// with the filter of `rhs`, OR if the memberships of `rhs` contain at least one bit set to 1 in common + /// with the filter of `self`. #[inline] pub const fn test(self, rhs: Self) -> bool { // NOTE: since const ops is not stable, we have to convert `Group` into u32 // to use & operator in const context. (self.memberships.bits() & rhs.filter.bits()) != 0 - && (rhs.memberships.bits() & self.filter.bits()) != 0 + || (rhs.memberships.bits() & self.filter.bits()) != 0 } } diff --git a/src/sampling/ray_sampling.rs b/src/sampling/ray_sampling.rs index 91fb10f..8e908d3 100644 --- a/src/sampling/ray_sampling.rs +++ b/src/sampling/ray_sampling.rs @@ -1,4 +1,4 @@ -use crate::math::{Isometry, Point, Real, Vector, DIM}; +use crate::math::{Real, Vector, DIM}; use parry::bounding_volume::{Aabb, BoundingVolume}; use parry::query::{Ray, RayCast}; @@ -9,8 +9,8 @@ use std::collections::HashSet; pub fn shape_surface_ray_sample( shape: &S, particle_rad: Real, -) -> Option>> { - let aabb = shape.compute_aabb(&Isometry::identity()); +) -> Option>> { + let aabb = shape.compute_aabb(&parry::math::Pose::IDENTITY); Some(surface_ray_sample(shape, &aabb, particle_rad)) } @@ -18,8 +18,8 @@ pub fn shape_surface_ray_sample( pub fn shape_volume_ray_sample( shape: &S, particle_rad: Real, -) -> Option>> { - let aabb = shape.compute_aabb(&Isometry::identity()); +) -> Option>> { + let aabb = shape.compute_aabb(&parry::math::Pose::IDENTITY); Some(volume_ray_sample(shape, &aabb, particle_rad)) } @@ -28,23 +28,25 @@ pub fn surface_ray_sample( shape: &S, volume: &Aabb, particle_rad: Real, -) -> Vec> { +) -> Vec> { let mut quantized_points = HashSet::new(); let subdivision_size = particle_rad * na::convert::<_, Real>(2.0); let volume = volume.loosened(subdivision_size); - let maxs = volume.maxs; - let origin = volume.mins + Vector::repeat(subdivision_size / na::convert::<_, Real>(2.0)); + let maxs: Vector = volume.maxs.into(); + let origin: Vector = (volume.mins + + parry::math::Vector::splat(subdivision_size / na::convert::<_, Real>(2.0))) + .into(); let mut curr = origin; - let mut perform_cast = |i, curr| { + let mut perform_cast = |i, curr: Vector| { let mut dir = Vector::zeros(); dir[i] = na::one::(); - let mut ray = Ray::new(curr, dir); + let mut ray = Ray::new(curr.into(), dir.into()); let mut entry_point = true; while let Some(toi) = shape.cast_local_ray(&ray, Real::MAX, false) { - let impact = ray.point_at(toi); + let impact: Vector = ray.point_at(toi).into(); let quantized_pt = quantize_point(&origin, &impact, subdivision_size, entry_point, i); let _ = quantized_points.insert(quantized_pt); ray.origin[i] += toi + subdivision_size / na::convert::<_, Real>(10.0); @@ -92,18 +94,20 @@ pub fn volume_ray_sample( shape: &S, volume: &Aabb, particle_rad: Real, -) -> Vec> { +) -> Vec> { let mut quantized_points = HashSet::new(); let subdivision_size = particle_rad * na::convert::<_, Real>(2.0); let volume = volume.loosened(subdivision_size); - let maxs = volume.maxs; - let origin = volume.mins + Vector::repeat(subdivision_size / na::convert::<_, Real>(2.0)); + let maxs: Vector = volume.maxs.into(); + let origin: Vector = (volume.mins + + parry::math::Vector::splat(subdivision_size / na::convert::<_, Real>(2.0))) + .into(); - let mut perform_cast = |i, curr| { + let mut perform_cast = |i, curr: Vector| { let mut dir = Vector::zeros(); dir[i] = na::one::(); - let mut ray = Ray::new(curr, dir); + let mut ray = Ray::new(curr.into(), dir.into()); let mut prev_impact = None; while let Some(toi) = shape.cast_local_ray(&ray, Real::MAX, false) { @@ -164,13 +168,13 @@ pub fn volume_ray_sample( } fn sample_segment( - origin: &Point, - start: &Point, + origin: &Vector, + start: &Vector, a: Real, b: Real, subdivision_size: Real, dimension: usize, - out: &mut HashSet>, + out: &mut HashSet>, ) { let mut quantized_pt = (start - origin).map(|e| { na::try_convert::<_, f64>(e / subdivision_size) @@ -191,28 +195,23 @@ fn sample_segment( } fn unquantize_points( - origin: &Point, + origin: &Vector, subdivision_size: Real, - quantized_points: &HashSet>, -) -> Vec> { + quantized_points: &HashSet>, +) -> Vec> { quantized_points .iter() - .map(|qpt| { - origin - + qpt - .coords - .map(|e| na::convert::<_, Real>(e as f64) * subdivision_size) - }) + .map(|qpt| origin + qpt.map(|e| na::convert::<_, Real>(e as f64) * subdivision_size)) .collect() } fn quantize_point( - origin: &Point, - point: &Point, + origin: &Vector, + point: &Vector, subdivision_size: Real, entry_point: bool, leading_dimension: usize, -) -> Point { +) -> Vector { let mut dpt = point - origin; for i in 0..DIM { if i == leading_dimension { diff --git a/src/solver/elasticity/becker2009_elasticity.rs b/src/solver/elasticity/becker2009_elasticity.rs index edd2953..b89af44 100644 --- a/src/solver/elasticity/becker2009_elasticity.rs +++ b/src/solver/elasticity/becker2009_elasticity.rs @@ -7,7 +7,7 @@ use approx::AbsDiffEq; use crate::geometry::{self, ParticlesContacts}; use crate::kernel::{CubicSplineKernel, Kernel}; -use crate::math::{Matrix, Point, Real, RotationMatrix, SpatialVector, Vector}; +use crate::math::{Matrix, Real, RotationMatrix, SpatialVector, Vector}; use crate::object::{Boundary, Fluid}; use crate::solver::NonPressureForce; use crate::TimestepManager; @@ -47,7 +47,7 @@ pub struct Becker2009Elasticity< d2: Real, nonlinear_strain: bool, volumes0: Vec, - positions0: Vec>, + positions0: Vec>, contacts0: ParticlesContacts, rotations: Vec>, deformation_gradient_tr: Vec>, diff --git a/src/solver/pressure/dfsph_solver.rs b/src/solver/pressure/dfsph_solver.rs index 74bcc8b..5b02663 100644 --- a/src/solver/pressure/dfsph_solver.rs +++ b/src/solver/pressure/dfsph_solver.rs @@ -13,6 +13,36 @@ use crate::object::{Boundary, Fluid}; use crate::solver::{helper, PressureSolver}; use crate::TimestepManager; +/// Tuning parameters for the DFSPH solver. +#[derive(Copy, Clone, Debug, PartialEq)] +pub struct DfsphParameters { + /// Minimum number of pressure iterations. + pub min_pressure_iter: usize, + /// Maximum number of pressure iterations. + pub max_pressure_iter: usize, + /// Maximum acceptable density error. + pub max_density_error: Real, + /// Minimum number of divergence iterations. + pub min_divergence_iter: usize, + /// Maximum number of divergence iterations. + pub max_divergence_iter: usize, + /// Maximum acceptable divergence error. + pub max_divergence_error: Real, +} + +impl Default for DfsphParameters { + fn default() -> Self { + Self { + min_pressure_iter: 1, + max_pressure_iter: 50, + max_density_error: na::convert::<_, Real>(0.05), + min_divergence_iter: 1, + max_divergence_iter: 50, + max_divergence_error: na::convert::<_, Real>(0.1), + } + } +} + /// A DFSPH (Divergence Free Smoothed Particle Hydrodynamics) pressure solver. pub struct DFSPHSolver< KernelDensity: Kernel = CubicSplineKernel, @@ -52,13 +82,14 @@ where { /// Initialize a new DFSPH pressure solver. pub fn new() -> Self { + let parameters = DfsphParameters::default(); Self { - min_pressure_iter: 1, - max_pressure_iter: 50, - max_density_error: na::convert::<_, Real>(0.05), - min_divergence_iter: 1, - max_divergence_iter: 50, - max_divergence_error: na::convert::<_, Real>(0.1), + min_pressure_iter: parameters.min_pressure_iter, + max_pressure_iter: parameters.max_pressure_iter, + max_density_error: parameters.max_density_error, + min_divergence_iter: parameters.min_divergence_iter, + max_divergence_iter: parameters.max_divergence_iter, + max_divergence_error: parameters.max_divergence_error, min_neighbors_for_divergence_solve: if DIM == 2 { 6 } else { 20 }, alphas: Vec::new(), densities: Vec::new(), @@ -69,6 +100,35 @@ where } } + /// Initialize a new DFSPH pressure solver with custom parameters. + pub fn with_parameters(parameters: DfsphParameters) -> Self { + let mut solver = Self::new(); + solver.set_parameters(parameters); + solver + } + + /// Current tuning parameters. + pub fn parameters(&self) -> DfsphParameters { + DfsphParameters { + min_pressure_iter: self.min_pressure_iter, + max_pressure_iter: self.max_pressure_iter, + max_density_error: self.max_density_error, + min_divergence_iter: self.min_divergence_iter, + max_divergence_iter: self.max_divergence_iter, + max_divergence_error: self.max_divergence_error, + } + } + + /// Set tuning parameters. + pub fn set_parameters(&mut self, parameters: DfsphParameters) { + self.max_pressure_iter = parameters.max_pressure_iter.max(1); + self.min_pressure_iter = parameters.min_pressure_iter.min(self.max_pressure_iter); + self.max_density_error = parameters.max_density_error; + self.max_divergence_iter = parameters.max_divergence_iter.max(1); + self.min_divergence_iter = parameters.min_divergence_iter.min(self.max_divergence_iter); + self.max_divergence_error = parameters.max_divergence_error; + } + fn compute_boundary_volumes( &mut self, boundary_boundary_contacts: &[ParticlesContacts], @@ -89,8 +149,11 @@ where denominator += c.weight; } - assert!(!denominator.is_zero()); - *volume = na::one::() / denominator; + if denominator.is_zero() { + *volume = na::zero(); + } else { + *volume = na::one::() / denominator; + } }) } } @@ -523,6 +586,15 @@ where KernelDensity: Kernel, KernelGradient: Kernel, { + fn dfsph_parameters(&self) -> Option { + Some(self.parameters()) + } + + fn set_dfsph_parameters(&mut self, parameters: DfsphParameters) -> bool { + self.set_parameters(parameters); + true + } + fn init_with_fluids(&mut self, fluids: &[Fluid]) { // Resize every buffer. self.alphas.resize(fluids.len(), Vec::new()); diff --git a/src/solver/pressure/iisph_solver.rs b/src/solver/pressure/iisph_solver.rs index 135cf13..9e3c8e5 100644 --- a/src/solver/pressure/iisph_solver.rs +++ b/src/solver/pressure/iisph_solver.rs @@ -83,8 +83,11 @@ where denominator += c.weight; } - assert!(!denominator.is_zero()); - *volume = na::one::() / denominator; + if denominator.is_zero() { + *volume = na::zero(); + } else { + *volume = na::one::() / denominator; + } }) } } diff --git a/src/solver/pressure/mod.rs b/src/solver/pressure/mod.rs index ecdf262..e419cbd 100644 --- a/src/solver/pressure/mod.rs +++ b/src/solver/pressure/mod.rs @@ -1,4 +1,4 @@ -pub use self::dfsph_solver::DFSPHSolver; +pub use self::dfsph_solver::{DFSPHSolver, DfsphParameters}; pub use self::iisph_solver::IISPHSolver; pub use self::pressure_solver::PressureSolver; diff --git a/src/solver/pressure/pressure_solver.rs b/src/solver/pressure/pressure_solver.rs index 2e47cac..ce4d8cb 100644 --- a/src/solver/pressure/pressure_solver.rs +++ b/src/solver/pressure/pressure_solver.rs @@ -2,10 +2,21 @@ use crate::counters::Counters; use crate::geometry::ContactManager; use crate::math::{Real, Vector}; use crate::object::{Boundary, Fluid}; +use crate::solver::DfsphParameters; use crate::TimestepManager; /// Trait implemented by pressure solvers. pub trait PressureSolver { + /// Current DFSPH tuning parameters, if this solver is DFSPH. + fn dfsph_parameters(&self) -> Option { + None + } + + /// Set DFSPH tuning parameters if this solver is DFSPH. + fn set_dfsph_parameters(&mut self, _parameters: DfsphParameters) -> bool { + false + } + /// Initialize this solver with the given fluids. fn init_with_fluids(&mut self, fluids: &[Fluid]); diff --git a/src/z_order.rs b/src/z_order.rs index 2b4cb51..de99561 100644 --- a/src/z_order.rs +++ b/src/z_order.rs @@ -1,4 +1,4 @@ -use crate::math::{Point, Real}; +use crate::math::{Real, Vector}; use num_traits::float::FloatCore; use std::cmp::Ordering; @@ -6,11 +6,10 @@ pub fn apply_permutation(permutation: &[usize], data: &[T]) -> Vec permutation.iter().map(|i| data[*i].clone()).collect() } -pub fn compute_points_z_order(points: &[Point]) -> Vec { +pub fn compute_points_z_order(points: &[Vector]) -> Vec { let mut indices: Vec<_> = (0..points.len()).collect(); indices.sort_unstable_by(|i, j| { - z_order_floats(points[*i].coords.as_slice(), points[*j].coords.as_slice()) - .unwrap_or(Ordering::Equal) + z_order_floats(points[*i].as_slice(), points[*j].as_slice()).unwrap_or(Ordering::Equal) }); indices }