diff --git a/Cargo.lock b/Cargo.lock
index 473bcf3..b23b141 100644
--- a/Cargo.lock
+++ b/Cargo.lock
@@ -136,7 +136,7 @@ dependencies = [
  "ndk-context",
  "ndk-sys 0.6.0+11769913",
  "num_enum",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -1022,6 +1022,19 @@ dependencies = [
  "syn",
 ]
 
+[[package]]
+name = "bevy_rich_text3d"
+version = "0.2.0"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "57601b9430ef61d3f17d74f3da0731460bf00fc465da4c45626214b8be1152f6"
+dependencies = [
+ "bevy",
+ "cosmic-text",
+ "sys-locale",
+ "thiserror 2.0.12",
+ "zeno 0.3.2",
+]
+
 [[package]]
 name = "bevy_scene"
 version = "0.15.3"
@@ -1235,6 +1248,7 @@ dependencies = [
  "bevy",
  "bevy_editor_cam",
  "bevy_rapier3d",
+ "bevy_rich_text3d",
  "bytemuck",
  "futures",
  "nalgebra",
@@ -1465,7 +1479,7 @@ dependencies = [
  "polling",
  "rustix",
  "slab",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -1965,7 +1979,7 @@ dependencies = [
  "encase_derive",
  "glam",
  "nalgebra",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -2459,7 +2473,7 @@ checksum = "c151a2a5ef800297b4e79efa4f4bec035c5f51d5ae587287c9b952bdf734cacd"
 dependencies = [
  "log",
  "presser",
- "thiserror",
+ "thiserror 1.0.69",
  "windows 0.58.0",
 ]
 
@@ -2682,7 +2696,7 @@ dependencies = [
  "combine",
  "jni-sys",
  "log",
- "thiserror",
+ "thiserror 1.0.69",
  "walkdir",
  "windows-sys 0.45.0",
 ]
@@ -2959,7 +2973,7 @@ dependencies = [
  "rustc-hash 1.1.0",
  "spirv",
  "termcolor",
- "thiserror",
+ "thiserror 1.0.69",
  "unicode-xid",
 ]
 
@@ -2978,7 +2992,7 @@ dependencies = [
  "regex",
  "regex-syntax 0.8.5",
  "rustc-hash 1.1.0",
- "thiserror",
+ "thiserror 1.0.69",
  "tracing",
  "unicode-ident",
 ]
@@ -3023,7 +3037,7 @@ dependencies = [
  "log",
  "ndk-sys 0.5.0+25.2.9519653",
  "num_enum",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -3038,7 +3052,7 @@ dependencies = [
  "ndk-sys 0.6.0+11769913",
  "num_enum",
  "raw-window-handle",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -3559,7 +3573,7 @@ dependencies = [
  "slab",
  "smallvec",
  "spade",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -3822,7 +3836,7 @@ dependencies = [
  "profiling",
  "rustc-hash 2.1.1",
  "simba",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -3955,7 +3969,7 @@ checksum = "6006a627c1a38d37f3d3a85c6575418cfe34a5392d60a686d0071e1c8d427acb"
 dependencies = [
  "cpal",
  "lewton",
- "thiserror",
+ "thiserror 1.0.69",
 ]
 
 [[package]]
@@ -4254,7 +4268,7 @@ checksum = "cbd59f3f359ddd2c95af4758c18270eddd9c730dde98598023cdabff472c2ca2"
 dependencies = [
  "skrifa",
  "yazi",
- "zeno",
+ "zeno 0.2.3",
 ]
 
 [[package]]
@@ -4318,7 +4332,16 @@ version = "1.0.69"
 source = "registry+https://github.com/rust-lang/crates.io-index"
 checksum = "b6aaf5339b578ea85b50e080feb250a3e8ae8cfcdff9a461c9ec2904bc923f52"
 dependencies = [
- "thiserror-impl",
+ "thiserror-impl 1.0.69",
+]
+
+[[package]]
+name = "thiserror"
+version = "2.0.12"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "567b8a2dae586314f7be2a752ec7474332959c6460e02bde30d702a66d488708"
+dependencies = [
+ "thiserror-impl 2.0.12",
 ]
 
 [[package]]
@@ -4332,6 +4355,17 @@ dependencies = [
  "syn",
 ]
 
+[[package]]
+name = "thiserror-impl"
+version = "2.0.12"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "7f7cf42b4507d8ea322120659672cf1b9dbb93f8f2d4ecfd6e51350ff5b17a1d"
+dependencies = [
+ "proc-macro2",
+ "quote",
+ "syn",
+]
+
 [[package]]
 name = "thread_local"
 version = "1.1.8"
@@ -4801,7 +4835,7 @@ dependencies = [
  "raw-window-handle",
  "rustc-hash 1.1.0",
  "smallvec",
- "thiserror",
+ "thiserror 1.0.69",
  "wgpu-hal",
  "wgpu-types",
 ]
@@ -4843,7 +4877,7 @@ dependencies = [
  "renderdoc-sys",
  "rustc-hash 1.1.0",
  "smallvec",
- "thiserror",
+ "thiserror 1.0.69",
  "wasm-bindgen",
  "web-sys",
  "wgpu-types",
@@ -5536,6 +5570,12 @@ version = "0.2.3"
 source = "registry+https://github.com/rust-lang/crates.io-index"
 checksum = "dd15f8e0dbb966fd9245e7498c7e9e5055d9e5c8b676b95bd67091cd11a1e697"
 
+[[package]]
+name = "zeno"
+version = "0.3.2"
+source = "registry+https://github.com/rust-lang/crates.io-index"
+checksum = "cc0de2315dc13d00e5df3cd6b8d2124a6eaec6a2d4b6a1c5f37b7efad17fcc17"
+
 [[package]]
 name = "zerocopy"
 version = "0.7.35"
diff --git a/Cargo.toml b/Cargo.toml
index a47571b..c868e4d 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -22,6 +22,7 @@ bevy_rapier3d = "0.29"
 
 [dev-dependencies]
 bevy_editor_cam = "0.5"
+bevy_rich_text3d = "0.2"
 
 [patch.crates-io]
 # rapier3d = { path = "../rapier/crates/rapier3d" }
diff --git a/examples/configurations.rs b/examples/configurations.rs
new file mode 100644
index 0000000..96d0d5b
--- /dev/null
+++ b/examples/configurations.rs
@@ -0,0 +1,364 @@
+use bevy::input::common_conditions;
+use bevy::prelude::*;
+use bevy::render::renderer::RenderDevice;
+use bevy_editor_cam::DefaultEditorCamPlugins;
+use bevy_editor_cam::prelude::EditorCam;
+use bevy_rapier3d::geometry::RapierColliderHandle;
+use bevy_rapier3d::plugin::ReadRapierContext;
+use bevy_rapier3d::prelude::{Collider, RigidBody};
+use bevy_rapier3d::render::RapierDebugRenderPlugin;
+use bevy_rich_text3d::{Text3d, Text3dBounds, Text3dPlugin, Text3dStyling, TextAtlas};
+use bevy_wgsparkl::components::MpmCouplingEnabled;
+use bevy_wgsparkl::resources::{AppState, PhysicsContext};
+use nalgebra::{Vector3, vector};
+use wgrapier3d::dynamics::body::{BodyCoupling, BodyCouplingEntry};
+use wgsparkl3d::models::DruckerPrager;
+use wgsparkl3d::solver::ParticlePhase;
+use wgsparkl3d::{
+    models::ElasticCoefficients,
+    pipeline::MpmData,
+    solver::{Particle, ParticleMassProps, SimulationParams},
+};
+
+pub fn main() {
+    App::new()
+        .add_plugins((DefaultPlugins, DefaultEditorCamPlugins))
+        .add_plugins(bevy_rapier3d::plugin::RapierPhysicsPlugin::<()>::default())
+        .add_plugins(RapierDebugRenderPlugin::default())
+        .add_plugins(bevy_wgsparkl::WgSparklPlugin)
+        .add_plugins(Text3dPlugin {
+            load_system_fonts: true,
+            ..Default::default()
+        })
+        .add_systems(PostUpdate, setup_mpm_particles)
+        .add_systems(
+            Update,
+            reset_scene.run_if(common_conditions::input_just_pressed(KeyCode::KeyR)),
+        )
+        .add_systems(Startup, setup_scene)
+        .run();
+}
+pub fn setup_scene(mut commands: Commands) {
+    commands.spawn((
+        Camera3d::default(),
+        EditorCam {
+            last_anchor_depth: 110f64,
+            ..Default::default()
+        },
+        Transform::from_xyz(0.0, 30.0, 100.0).looking_at(Vec3::new(0.0, 10.0, 0.0), Vec3::Y),
+    ));
+    /*
+     * Ground
+     */
+    let ground_size = 1200.1;
+    let ground_height = 2.0;
+
+    commands.spawn((
+        Transform::from_xyz(0.0, -ground_height, 0.0),
+        Collider::cuboid(ground_size, ground_height, ground_size),
+        RigidBody::Fixed,
+        MpmCouplingEnabled,
+    ));
+}
+
+pub fn reset_scene(mut commands: Commands, mut app_state: ResMut<AppState>) {
+    app_state.restarting = true;
+    app_state.particles_initialized = false;
+    commands.remove_resource::<PhysicsContext>();
+}
+
+pub fn setup_mpm_particles(
+    mut commands: Commands,
+    device: Res<RenderDevice>,
+    mut app_state: ResMut<AppState>,
+    rapier: ReadRapierContext,
+    coupling: Query<&RapierColliderHandle, With<MpmCouplingEnabled>>,
+    mut standard_materials: ResMut<Assets<StandardMaterial>>,
+) {
+    if rapier.rapier_context.get_single().is_err() {
+        return; // Rapier isn’t initialized yet.
+    }
+
+    let rapier = rapier.single();
+
+    if rapier.colliders.colliders.is_empty() {
+        return; // Rapier isn’t initialized yet.
+    }
+
+    if app_state.particles_initialized {
+        return; // Already initialized.
+    }
+
+    let grid_size_x = 10;
+    let grid_size_y = 10;
+    let grid_size_z = 10;
+    let num_particles = grid_size_x * grid_size_y * grid_size_z;
+
+    let particle_positions = (0..num_particles)
+        .map(|i| {
+            let x = i % grid_size_x;
+            let y = (i / grid_size_x) % grid_size_y;
+            let z = (i / (grid_size_x * grid_size_y)) % grid_size_z;
+            Vector3::new(x as f32, y as f32 + 1f32, z as f32)
+        })
+        .collect::<Vec<_>>();
+
+    app_state.particles_initialized = true;
+
+    let coupling: Vec<_> = coupling
+        .iter()
+        .map(|co_handle| {
+            let co = &rapier.colliders.colliders[co_handle.0];
+            println!("Coupled collider: {:?}", co.shape().shape_type());
+            println!(
+                "Coupled collider pose: {:?}",
+                co.position().translation.vector
+            );
+            let rb_handle = co.parent().unwrap();
+            BodyCouplingEntry {
+                body: rb_handle,
+                collider: co_handle.0,
+                mode: BodyCoupling::OneWay, // TODO: try out two-ways for the particles to affect the rigid bodies.
+            }
+        })
+        .collect();
+
+    let device = device.wgpu_device();
+
+    if !app_state.restarting {
+        app_state.num_substeps = 16;
+        app_state.gravity_factor = 1.0;
+    };
+
+    // Text material.
+    let mat = standard_materials.add(StandardMaterial {
+        base_color_texture: Some(TextAtlas::DEFAULT_IMAGE.clone_weak()),
+        alpha_mode: AlphaMode::Blend,
+        unlit: true,
+        ..Default::default()
+    });
+
+    let params = SimulationParams {
+        gravity: vector![0.0, -9.81, 0.0] * app_state.gravity_factor,
+        dt: (1.0 / 60.0) / (app_state.num_substeps as f32),
+    };
+
+    let cell_width = 1f32;
+    let mut particles = vec![];
+    let mut configurations = vec![];
+    let get_position_for_line = |z: f32| -> bevy::math::Vec3 {
+        bevy::math::Vec3::new(
+            -2f32 * grid_size_x as f32 * 3f32 * 2f32,
+            5f32,
+            z * grid_size_z as f32 * 0.7f32 * 3f32 * 2f32 + grid_size_z as f32 / 2f32,
+        )
+    };
+    let mut display_text_at_world_pos = |world_pos: bevy::math::Vec3, text: String| {
+        commands.spawn((
+            Text3d::new(text),
+            Text3dStyling {
+                size: 40.,
+                color: Srgba::new(1., 1., 1., 1.),
+                ..Default::default()
+            },
+            Text3dBounds { width: 500. },
+            Mesh3d::default(),
+            MeshMaterial3d(mat.clone()),
+            Transform::from_translation(world_pos)
+                .with_scale(Vec3::splat(0.1))
+                .with_rotation(Quat::from_rotation_x(-90f32.to_radians())),
+        ));
+    };
+    {
+        let young_modulus = 1_000_000_000.0;
+        let z = -1f32;
+        display_text_at_world_pos(
+            get_position_for_line(z),
+            format!(
+                "With plasticity.\nmodulus = {}M",
+                young_modulus / 1_000_000.0
+            ),
+        );
+        // line with plasticity, varying poisson
+        for x in -1..2 {
+            let poisson_ratio = match x {
+                -1 => 0.0,
+                0 => 0.2,
+                1 => 0.4,
+                _ => unreachable!(),
+            };
+            let model = ElasticCoefficients::from_young_modulus(young_modulus, poisson_ratio);
+            let plasticity = Some(DruckerPrager {
+                h0: 35.0f32.to_radians(),
+                h1: 9.0f32.to_radians(),
+                h2: 0.2,
+                h3: 10.0f32.to_radians(),
+                ..DruckerPrager::new(model.lambda, model.mu)
+            });
+            configurations.push(ParticlesConfiguration {
+                coords: IVec2::new(x, z as i32),
+                density: 3700f32,
+                model,
+                plasticity,
+                phase: None,
+                description: format!("poisson: {}", poisson_ratio),
+            });
+        }
+    }
+
+    {
+        let poisson_ratio = 0f32;
+        let z = 0f32;
+        display_text_at_world_pos(
+            get_position_for_line(z),
+            format!("With plasticity.\npoisson = {}", poisson_ratio),
+        );
+        // line with plasticity, varying young modulus
+        for x in -1..2 {
+            let young_modulus = match x {
+                -1 => 1_000_000.0,
+                0 => 10_000_000.0,
+                1 => 100_000_000.0,
+                _ => unreachable!(),
+            };
+            let model = ElasticCoefficients::from_young_modulus(young_modulus, poisson_ratio);
+            let plasticity = Some(DruckerPrager {
+                h0: 35.0f32.to_radians(),
+                h1: 9.0f32.to_radians(),
+                h2: 0.2,
+                h3: 10.0f32.to_radians(),
+                ..DruckerPrager::new(model.lambda, model.mu)
+            });
+            configurations.push(ParticlesConfiguration {
+                coords: IVec2::new(x, z as i32),
+                density: 3700f32,
+                model,
+                plasticity,
+                phase: None,
+                description: format!("modulus: {}M", young_modulus / 1_000_000f32),
+            });
+        }
+    }
+
+    {
+        let poisson_ratio = 0f32;
+        let z = 1f32;
+        display_text_at_world_pos(
+            get_position_for_line(z),
+            format!("Without plasticity.\npoisson = {}", poisson_ratio),
+        );
+        // line without plasticity, varying young modulus
+        for x in -1..2 {
+            let young_modulus = match x {
+                -1 => 1_000_000.0,
+                0 => 50_000_000.0,
+                1 => 200_000_000.0,
+                _ => unreachable!(),
+            };
+            let model = ElasticCoefficients::from_young_modulus(young_modulus, poisson_ratio);
+            configurations.push(ParticlesConfiguration {
+                coords: IVec2::new(x, z as i32),
+                density: 3700f32,
+                model,
+                plasticity: None,
+                phase: Some(ParticlePhase {
+                    phase: 1.0,
+                    max_stretch: f32::MAX,
+                }),
+                description: format!("modulus: {}M", young_modulus / 1_000_000.0),
+            });
+        }
+    }
+    {
+        let young_modulus = 1_000_000.0;
+        let z = 2f32;
+        display_text_at_world_pos(
+            get_position_for_line(z),
+            format!(
+                "Without plasticity.\nmodulus = {}M",
+                young_modulus / 1_000_000.0
+            ),
+        );
+        // line without plasticity, varying poisson_ratio
+        for x in -1..2 {
+            let poisson_ratio = match x {
+                -1 => -0.2f32,
+                0 => 0.3,
+                1 => 0.48,
+                _ => unreachable!(),
+            };
+            let model = ElasticCoefficients::from_young_modulus(young_modulus, poisson_ratio);
+            configurations.push(ParticlesConfiguration {
+                coords: IVec2::new(x, z as i32),
+                density: 3700f32,
+                model,
+                plasticity: None,
+                phase: Some(ParticlePhase {
+                    phase: 1.0,
+                    max_stretch: f32::MAX,
+                }),
+                description: format!("poisson: {}", poisson_ratio),
+            });
+        }
+    }
+
+    for c in configurations.iter() {
+        let x = c.coords.x as f32 * 3f32;
+        let z = c.coords.y as f32 * 3f32;
+        let offset = vector![
+            x * grid_size_x as f32,
+            3f32,
+            z * grid_size_z as f32 * 0.7f32
+        ] * 2f32;
+        for particle in &particle_positions {
+            let position = vector![particle.x, particle.y, particle.z];
+
+            let particle_size = vector![1.0, 1.0, 1.0];
+            let volume = particle_size.x * particle_size.y * particle_size.z;
+            let density = c.density;
+            particles.push(Particle {
+                position: nalgebra::Rotation::from_axis_angle(
+                    &Vector3::z_axis(),
+                    1f32.to_radians(),
+                ) * vector![position.x, position.y, position.z]
+                    + offset,
+                velocity: Vector3::zeros(),
+                volume: ParticleMassProps::new(density * volume, volume.cbrt() / 2.0),
+                model: c.model,
+                plasticity: c.plasticity,
+                phase: c.phase,
+            });
+        }
+
+        display_text_at_world_pos(
+            Vec3::new(offset.x, 5f32, offset.z + 10f32 + grid_size_z as f32),
+            c.description.clone(),
+        );
+    }
+
+    println!("Number of simulated particles: {}", particles.len());
+
+    println!("Coupled: {}", coupling.len());
+
+    let data = MpmData::with_select_coupling(
+        device,
+        params,
+        &particles,
+        &rapier.rigidbody_set.bodies,
+        &rapier.colliders.colliders,
+        coupling,
+        cell_width,
+        60_000,
+    );
+    commands.insert_resource(PhysicsContext { data, particles });
+}
+
+#[derive(Debug)]
+pub struct ParticlesConfiguration {
+    pub coords: IVec2,
+    pub density: f32,
+    pub model: ElasticCoefficients,
+    pub plasticity: Option<DruckerPrager>,
+    pub phase: Option<ParticlePhase>,
+    pub description: String,
+}