diff --git a/12_pbd_cloth/constraints/bending.py b/12_pbd_cloth/constraints/bending.py index 6417ddb..2e65cab 100644 --- a/12_pbd_cloth/constraints/bending.py +++ b/12_pbd_cloth/constraints/bending.py @@ -5,7 +5,7 @@ def solve_bending_constraints( compliance: ti.f64, dt: ti.f64, num_bending_constraints: ti.i32, pos: ti.template(), bending_ids: ti.template(), bending_lengths: ti.template(), - inv_mass: ti.template() + inv_mass: ti.template(), lambdas: ti.template() ): alpha = compliance / (dt * dt) for i in range(num_bending_constraints): @@ -19,6 +19,7 @@ def solve_bending_constraints( if dist == 0.0: continue grad = delta / dist C = dist - bending_lengths[i] - s = -C / (w_sum + alpha) - pos[id0] += s * w0 * grad - pos[id1] -= s * w1 * grad + dlambda = -(C + alpha * lambdas[i]) / (w_sum + alpha) + lambdas[i] += dlambda + pos[id0] += dlambda * w0 * grad + pos[id1] -= dlambda * w1 * grad diff --git a/12_pbd_cloth/constraints/stretching.py b/12_pbd_cloth/constraints/stretching.py index ebdce96..4c00eb8 100644 --- a/12_pbd_cloth/constraints/stretching.py +++ b/12_pbd_cloth/constraints/stretching.py @@ -5,7 +5,7 @@ def solve_stretching_constraints( compliance: ti.f64, dt: ti.f64, num_stretching_constraints: ti.i32, pos: ti.template(), stretching_ids: ti.template(), stretching_lengths: ti.template(), - inv_mass: ti.template() + inv_mass: ti.template(), lambdas: ti.template() ): alpha = compliance / (dt * dt) for i in range(num_stretching_constraints): @@ -19,6 +19,7 @@ def solve_stretching_constraints( if dist == 0.0: continue grad = delta / dist C = dist - stretching_lengths[i] - s = -C / (w_sum + alpha) - pos[id0] += s * w0 * grad - pos[id1] -= s * w1 * grad + dlambda = -(C + alpha * lambdas[i]) / (w_sum + alpha) + lambdas[i] += dlambda + pos[id0] += dlambda * w0 * grad + pos[id1] -= dlambda * w1 * grad diff --git a/12_pbd_cloth/imgui.ini b/12_pbd_cloth/imgui.ini new file mode 100644 index 0000000..ea675fc --- /dev/null +++ b/12_pbd_cloth/imgui.ini @@ -0,0 +1,10 @@ +[Window][Debug##Default] +Pos=60,60 +Size=400,400 +Collapsed=0 + +[Window][Controls] +Pos=51,51 +Size=307,409 +Collapsed=0 + diff --git a/12_pbd_cloth/main.py b/12_pbd_cloth/main.py index a136491..6426228 100644 --- a/12_pbd_cloth/main.py +++ b/12_pbd_cloth/main.py @@ -78,7 +78,7 @@ def save(self): dt = 1.0 / 60.0 num_substeps = 15 sdt = dt / num_substeps -solver_iterations = 1 +solver_iterations = 1 # can increase export_enabled = False export_frame_count = 0 @@ -112,6 +112,10 @@ def save(self): stretching_lengths = ti.field(dtype=ti.f64, shape=num_stretching_constraints) bending_lengths = ti.field(dtype=ti.f64, shape=num_bending_constraints) +# XPBD accumulated lambdas (one per constraint, reset each substep) +stretching_lambdas = ti.field(dtype=ti.f64, shape=num_stretching_constraints) +bending_lambdas = ti.field(dtype=ti.f64, shape=num_bending_constraints) + # Visualization fields ground_vertices = ti.Vector.field(3, dtype=ti.f64, shape=4) ground_indices = ti.field(ti.i32, shape=6) @@ -133,12 +137,28 @@ def save(self): # Simulation Substep Function # ============================================================================ +@ti.kernel +def reset_constraint_lambdas( + num_s: ti.i32, num_b: ti.i32, + s_lambdas: ti.template(), b_lambdas: ti.template() +): + for i in range(num_s): + s_lambdas[i] = 0.0 + for i in range(num_b): + b_lambdas[i] = 0.0 + + def substep(grab_id, grab_x, grab_y, grab_z): """Perform one simulation substep.""" pre_solve(sdt, num_particles, gravity, pos, prev_pos, vel, inv_mass) + # XPBD: reset lambdas once per substep (they accumulate across solver iterations) + reset_constraint_lambdas( + num_stretching_constraints, num_bending_constraints, + stretching_lambdas, bending_lambdas + ) for _ in range(solver_iterations): - solve_stretching_constraints(stretching_compliance, sdt, num_stretching_constraints, pos, stretching_ids, stretching_lengths, inv_mass) - solve_bending_constraints(bending_compliance, sdt, num_bending_constraints, pos, bending_ids, bending_lengths, inv_mass) + solve_stretching_constraints(stretching_compliance, sdt, num_stretching_constraints, pos, stretching_ids, stretching_lengths, inv_mass, stretching_lambdas) + solve_bending_constraints(bending_compliance, sdt, num_bending_constraints, pos, bending_ids, bending_lengths, inv_mass, bending_lambdas) apply_grab(grab_id, grab_x, grab_y, grab_z, pos, vel, grab_indicator_pos) post_solve(sdt, num_particles, pos, prev_pos, vel, inv_mass) diff --git a/12_pbd_cloth/utils/cloth_utils.py b/12_pbd_cloth/utils/cloth_utils.py index 3f81fd8..dbea5e8 100644 --- a/12_pbd_cloth/utils/cloth_utils.py +++ b/12_pbd_cloth/utils/cloth_utils.py @@ -23,13 +23,13 @@ def initialize_colors( ): for i in range(num_particles): - x_coord = int((pos[i].x + 1.0) * 5) % 2 - z_coord = int((pos[i].z + 1.0) * 5) % 2 + x_coord = int((pos[i].x + 1.0) * 8) % 2 + y_coord = int((pos[i].y + 1.0) * 8) % 2 - if (x_coord + z_coord) % 2 == 0: - vertex_colors[i] = ti.Vector([1.0, 1.0, 1.0]) + if (x_coord + y_coord) % 2 == 0: + vertex_colors[i] = ti.Vector([0.15, 0.35, 0.75]) # deep blue else: - vertex_colors[i] = ti.Vector([0.1, 0.1, 0.1]) + vertex_colors[i] = ti.Vector([0.95, 0.92, 0.85]) # light cream @ti.kernel diff --git a/13_pbd_mesh/constraints/edge.py b/13_pbd_mesh/constraints/edge.py index ccf3228..43d8dd7 100644 --- a/13_pbd_mesh/constraints/edge.py +++ b/13_pbd_mesh/constraints/edge.py @@ -4,7 +4,8 @@ @ti.kernel def solve_edges( compliance: ti.f64, dt: ti.f64, num_edges: ti.i32, - pos: ti.template(), edge_ids: ti.template(), edge_lengths: ti.template(), inv_mass: ti.template() + pos: ti.template(), edge_ids: ti.template(), edge_lengths: ti.template(), + inv_mass: ti.template(), lambdas: ti.template() ): alpha = compliance / (dt * dt) for i in range(num_edges): @@ -17,6 +18,7 @@ def solve_edges( if dist == 0.0: continue grad = delta / dist C = dist - edge_lengths[i] - s = -C / (w_sum + alpha) - pos[id0] += s * w0 * grad - pos[id1] -= s * w1 * grad + dlambda = -(C + alpha * lambdas[i]) / (w_sum + alpha) + lambdas[i] += dlambda + pos[id0] += dlambda * w0 * grad + pos[id1] -= dlambda * w1 * grad diff --git a/13_pbd_mesh/constraints/volume.py b/13_pbd_mesh/constraints/volume.py index 493d596..74a4735 100644 --- a/13_pbd_mesh/constraints/volume.py +++ b/13_pbd_mesh/constraints/volume.py @@ -12,7 +12,7 @@ def solve_volumes( compliance: ti.f64, dt: ti.f64, num_tets: ti.i32, pos: ti.template(), tet_ids: ti.template(), rest_vol: ti.template(), inv_mass: ti.template(), - vol_id_order: ti.template() + vol_id_order: ti.template(), lambdas: ti.template() ): alpha = compliance / (dt * dt) for i in range(num_tets): @@ -27,6 +27,7 @@ def solve_volumes( w_sum += inv_mass[p_indices[j]] * grad.norm_sqr() if w_sum == 0.0: continue C = get_tet_volume(p_indices, pos) - rest_vol[i] - s = -C / (w_sum + alpha) + dlambda = -(C + alpha * lambdas[i]) / (w_sum + alpha) + lambdas[i] += dlambda for j in ti.static(range(4)): - pos[p_indices[j]] += s * inv_mass[p_indices[j]] * grads[j, :] + pos[p_indices[j]] += dlambda * inv_mass[p_indices[j]] * grads[j, :] diff --git a/13_pbd_mesh/imgui.ini b/13_pbd_mesh/imgui.ini new file mode 100644 index 0000000..b7cc025 --- /dev/null +++ b/13_pbd_mesh/imgui.ini @@ -0,0 +1,10 @@ +[Window][Debug##Default] +Pos=60,60 +Size=400,400 +Collapsed=0 + +[Window][Controls] +Pos=51,51 +Size=307,256 +Collapsed=0 + diff --git a/13_pbd_mesh/main.py b/13_pbd_mesh/main.py index 4d4369c..8af65a4 100644 --- a/13_pbd_mesh/main.py +++ b/13_pbd_mesh/main.py @@ -154,6 +154,10 @@ def save(self): tet_ids = ti.field(ti.i32, shape=(num_tets, 4)) edge_ids = ti.field(ti.i32, shape=(num_edges, 2)) +# XPBD accumulated lambdas (one per constraint, reset each substep) +edge_lambdas = ti.field(dtype=ti.f64, shape=num_edges) +vol_lambdas = ti.field(dtype=ti.f64, shape=num_tets) + # Visual mesh fields vis_mesh_rest_pos = ti.Vector.field(3, dtype=ti.f64, shape=num_vis_verts) vis_mesh_pos = ti.Vector.field(3, dtype=ti.f64, shape=num_vis_verts) @@ -194,11 +198,24 @@ def save(self): # Simulation Substep Function # ============================================================================ +@ti.kernel +def reset_constraint_lambdas( + num_e: ti.i32, num_v: ti.i32, + e_lambdas: ti.template(), v_lambdas: ti.template() +): + for i in range(num_e): + e_lambdas[i] = 0.0 + for i in range(num_v): + v_lambdas[i] = 0.0 + + def substep(): pre_solve(sdt, 1, num_particles, gravity, box_min, box_max, pos, prev_pos, vel, inv_mass) + # XPBD: reset lambdas once per substep (they accumulate across solver iterations) + reset_constraint_lambdas(num_edges, num_tets, edge_lambdas, vol_lambdas) for _ in range(solver_iterations): - solve_edges(edge_compliance, sdt, num_edges, pos, edge_ids, edge_lengths, inv_mass) - solve_volumes(vol_compliance, sdt, num_tets, pos, tet_ids, rest_vol, inv_mass, vol_id_order) + solve_edges(edge_compliance, sdt, num_edges, pos, edge_ids, edge_lengths, inv_mass, edge_lambdas) + solve_volumes(vol_compliance, sdt, num_tets, pos, tet_ids, rest_vol, inv_mass, vol_id_order, vol_lambdas) post_solve(sdt, num_particles, pos, prev_pos, vel, inv_mass) @@ -284,6 +301,7 @@ def substep(): scene.mesh(vis_mesh_pos, indices=vis_mesh_indices, per_vertex_color=vis_mesh_colors) scene.mesh(ground_vertices, indices=ground_indices, per_vertex_color=ground_colors) + canvas.set_background_color((1.0, 1.0, 1.0)) # white canvas.scene(scene) window.show()