refactor:rigid-body rod contact friction examples

examples/RigidBodyCases/RodRigidBodyContact/post_processing.py

@@ -728,8 +728,9 @@ def plot_video_with_surface(
 
 
 def plot_force_vs_energy(
-    force,
+    normalized_force,
     total_final_energy,
+    friction_coefficient,
     filename="energy_vs_force.png",
     SAVE_FIGURE=False,
 ):
@@ -739,12 +740,13 @@ def plot_force_vs_energy(
     axs.append(plt.subplot2grid((1, 1), (0, 0)))
 
     axs[0].plot(
-        force,
+        normalized_force,
         total_final_energy,
         linewidth=3,
     )
+    plt.axvline(x=friction_coefficient, linewidth=3, color="r", label="threshold")
     axs[0].set_ylabel("total energy", fontsize=20)
-    axs[0].set_xlabel("force ", fontsize=20)
+    axs[0].set_xlabel("normalized force", fontsize=20)
 
     plt.tight_layout()
     # fig.align_ylabels()

examples/RigidBodyCases/RodRigidBodyContact/rod_cylinder_contact_friction.py

@@ -6,7 +6,9 @@
 from post_processing import plot_velocity, plot_video_with_surface
 
 
-def rod_cylinder_contact_friction_case(pulling_force=0.0, POST_PROCESSING=False):
+def rod_cylinder_contact_friction_case(
+    force_coefficient=0.1, normal_force_mag=10, POST_PROCESSING=False
+):
     class RodCylinderParallelContact(
         BaseSystemCollection, Constraints, Connections, CallBacks, Forcing
     ):
@@ -15,7 +17,7 @@ class RodCylinderParallelContact(
     rod_cylinder_parallel_contact_simulator = RodCylinderParallelContact()
 
     # time step etc
-    final_time = 10.0
+    final_time = 20.0
     time_step = 1e-4
     total_steps = int(final_time / time_step) + 1
     rendering_fps = 30  # 20 * 1e1
@@ -47,18 +49,18 @@ class RodCylinderParallelContact(
     )
 
     rod_cylinder_parallel_contact_simulator.append(rod)
-    # rod.velocity_collection[0, :] -= 0.2
 
-    # Apply uniform forces on the rod
+    # Push the rod towards the cylinder to make sure contact is there
+    normal_force_direction = np.array([-1.0, 0.0, 0.0])
     rod_cylinder_parallel_contact_simulator.add_forcing_to(rod).using(
-        UniformForces, force=1.0 * pulling_force, direction=direction
+        UniformForces, force=normal_force_mag, direction=normal_force_direction
     )
-    # Push the rod towards the cylinder to make sure contact is there
+    # Apply uniform forces on the rod
     rod_cylinder_parallel_contact_simulator.add_forcing_to(rod).using(
-        UniformForces, force=0.1, direction=np.array([-1.0, 0.0, 0.0])
+        UniformForces, force=normal_force_mag * force_coefficient, direction=direction
     )
 
-    cylinder_height = 3 * base_length
+    cylinder_height = 8 * base_length
     cylinder_radius = base_radius
 
     cylinder_start = start + np.array([-1.0, 0.0, 0.0]) * 2 * base_radius
@@ -83,10 +85,10 @@ class RodCylinderParallelContact(
     # Add contact between rigid body and rod
     rod_cylinder_parallel_contact_simulator.connect(rod, rigid_body).using(
         ExternalContact,
-        k=1e4,
-        nu=10,
-        velocity_damping_coefficient=10,
-        kinetic_friction_coefficient=10,
+        k=1e5,
+        nu=100,
+        velocity_damping_coefficient=1e5,
+        friction_coefficient=0.5,
     )
 
     # Add callbacks
@@ -259,21 +261,3 @@ def make_callback(self, system, time, current_step: int):
     )
 
     return total_final_energy
-
-
-if __name__ == "__main__":
-    import multiprocessing as mp
-    from examples.RigidBodyCases.RodRigidBodyContact.post_processing import (
-        plot_force_vs_energy,
-    )
-
-    # total_energy = rod_cylinder_contact_friction_case(pulling_force=0.5, POST_PROCESSING=True)
-
-    force = list([(float(x)) / 100.0 for x in range(0, 90, 5)])
-
-    with mp.Pool(mp.cpu_count()) as pool:
-        result_total_energy = pool.map(rod_cylinder_contact_friction_case, force)
-
-    plot_force_vs_energy(
-        force, result_total_energy, filename="rod_energy_vs_force.png", SAVE_FIGURE=True
-    )

examples/RigidBodyCases/RodRigidBodyContact/rod_cylinder_contact_friction_case.py

@@ -0,0 +1,8 @@
+if __name__ == "__main__":
+    from examples.RigidBodyCases.RodRigidBodyContact.rod_cylinder_contact_friction import (
+        rod_cylinder_contact_friction_case,
+    )
+
+    total_energy = rod_cylinder_contact_friction_case(
+        force_coefficient=0.1, normal_force_mag=10, POST_PROCESSING=True
+    )

examples/RigidBodyCases/RodRigidBodyContact/rod_cylinder_contact_friction_phase_space.py

@@ -0,0 +1,25 @@
+if __name__ == "__main__":
+    import multiprocessing as mp
+    from examples.RigidBodyCases.RodRigidBodyContact.rod_cylinder_contact_friction import (
+        rod_cylinder_contact_friction_case,
+    )
+    from examples.RigidBodyCases.RodRigidBodyContact.post_processing import (
+        plot_force_vs_energy,
+    )
+
+    # total_energy = rod_cylinder_contact_friction_case(friction_coefficient=1.0, normal_force_mag=10, velocity_damping_coefficient=1E4, POST_PROCESSING=True)
+
+    friction_coefficient = list([(float(x)) / 100.0 for x in range(0, 100, 5)])
+
+    with mp.Pool(mp.cpu_count()) as pool:
+        result_total_energy = pool.map(
+            rod_cylinder_contact_friction_case, friction_coefficient
+        )
+
+    plot_force_vs_energy(
+        friction_coefficient,
+        result_total_energy,
+        friction_coefficient=0.5,
+        filename="rod_energy_vs_force.png",
+        SAVE_FIGURE=True,
+    )