Improve ContactResult comments and test

- Improve comments in elasticCollisionResponse.jl.
- Improve plot definitions in BouncingSphereContactResults.jl.

Author: GerhardHippmann

src/Composition/responseCalculation/elasticCollisionResponse.jl

@@ -138,50 +138,37 @@ function responseCalculation(material::ElasticContactPairResponseMaterial, obj1:
     vsmall   = F(material.vsmall)
     wsmall   = F(material.wsmall)
 
-    ### signed velocity and relative velocity ####
     # Contact points and distances to local part frame (in world frame)
-    # r_rel1 ... position vector from obj1 to average contact point, resolved in world
-    # r_rel2 ... position vector from obj2 to average contact point, resolved in world
-    r_rel1 = rContact - obj1.r_abs
-    r_rel2 = rContact - obj2.r_abs
+    r_rel1 = rContact - obj1.r_abs  # position vector from obj1 to average contact point, resolved in world
+    r_rel2 = rContact - obj2.r_abs  # position vector from obj2 to average contact point, resolved in world
 
     # Velocities and angular velocities of contact frames in world frame
-
-    w1 = obj1.R_abs'*obj1.w
-    w2 = obj2.R_abs'*obj2.w
-    v1 = obj1.v0 + cross(w1,r_rel1)
-    v2 = obj2.v0 + cross(w2,r_rel2)
+    w1 = obj1.R_abs'*obj1.w           # absolute angular velocity vector of obj1, resolved in world
+    w2 = obj2.R_abs'*obj2.w           # absolute angular velocity vector of obj2, resolved in world
+    v1 = obj1.v0 + cross(w1, r_rel1)  # absolute velocity vector of obj1 at average contact point, resolved in world
+    v2 = obj2.v0 + cross(w2, r_rel2)  # absolute velocity vector of obj2 at average contact point, resolved in world
 
     # Velocities and angular velocities in normal and tangential direction
-    # w_rel   ... angular velocity vector of obj2 relative to obj1, resolved in world
-    # e_w_reg ... regularised direction vector of w_rel, resolved in world
-    # v_rel   ... velocity vector at average contact point on obj2 relative to obj1, resolved in world
-    w_rel   = w2 - w1
-    e_w_reg = w_rel/Modia3D.regularize(norm(w_rel),wsmall)
-    v_rel   = v2 - v1
-
-    # e_n       ... points from contact point on obj2 to contact point on obj1, i.e. normal direction into obj2
-    # s         ... signed distance of the contact points, negative during contact
-    # delta_dot ... signed relative velocity in normal direction, positive for expansion (sign inconsistent with delta!)
-    # v_t       ... tangential velocity vector at average contact point on obj2 relative to obj1, resolved in world
-    # e_t_reg   ... regularised direction vector of v_t, resolved in world
-    # delta     ... penetration, positive during contact
-    delta_dot = dot(v_rel,e_n)
-    v_t       = v_rel - delta_dot*e_n
-    e_t_reg   = v_t/Modia3D.regularize(norm(v_t),vsmall)
-    delta     = -s
-
+    w_rel   = w2 - w1                                        # angular velocity vector of obj2 relative to obj1, resolved in world
+    e_w_reg = w_rel/Modia3D.regularize(norm(w_rel), wsmall)  # regularized direction vector of w_rel, resolved in world
+    v_rel   = v2 - v1                                        # velocity vector of obj2 relative to obj1 at average contact point, resolved in world
+
+    # e_n points from contact point on obj2 to contact point on obj1, i.e. normal direction into obj2
+    # s = signed distance of the contact points, negative during contact
+    delta_dot  = dot(v_rel, e_n)                            # signed relative velocity in normal direction, positive for expansion (sign inconsistent with delta!)
+    v_t        = v_rel - delta_dot*e_n                      # tangential velocity vector of obj2 relative to obj1 at average contact point, resolved in world
+    e_t_reg    = v_t/Modia3D.regularize(norm(v_t), vsmall)  # regularized direction vector of v_t, resolved in world
+    delta      = -s                                         # penetration, positive during contact
     delta_comp = delta * sqrt(abs(delta))
 
-    #fn = -max(F(0.0), c_res * c_geo * delta_comp * (1 - d_res*delta_dot) )
     fn = -c_res * c_geo * delta_comp * (1 - d_res*delta_dot)  # normal force, negative for pressure!
     ft = -mu_k*fn*e_t_reg                                     # tangential force vector acting on obj1, resolved in world
-    f1 = fn * e_n + ft                                        # total force vector acting on obj1, resolved in world
-    f2 = -f1                                                  # total force vector acting on obj2, resolved in world
+    f1 = fn * e_n + ft                                        # total force vector acting at average contact point on obj1, resolved in world
+    f2 = -f1                                                  # total force vector acting at average contact point on obj2, resolved in world
 
     tau = -mu_r * mu_r_geo * fn * e_w_reg                     # torque vector acting at average contact point on obj1, resolved in world
-    t1  = cross(r_rel1,f1) + tau                              # torque vector acting at obj1, resolved in world
-    t2  = cross(r_rel2,f2) - tau                              # torque vector acting at obj2, resolved in world
+    t1  = cross(r_rel1, f1) + tau                             # total torque vector acting at obj1, resolved in world
+    t2  = cross(r_rel2, f2) - tau                             # total torque vector acting at obj2, resolved in world
 
     return (f1, f2, t1, t2,
             ContactResults(delta,        # normal contact penetration (positive during contact)

test/Collision/BouncingSphereContactResults.jl

@@ -39,7 +39,8 @@ simulate!(bouncingSphere, stopTime=stopTime, tolerance=tolerance, dtmax=dtmax, l
           requiredFinalStates_atol=1e-6, requiredFinalStates=requiredFinalStates)
 
 @usingModiaPlot
-plot(bouncingSphere, ["result.penetration", "result.penetrationVelocity", "result.tangentialVelocity", "result.angularVelocity", "result.normalForce", "result.tangentialForce", "result.torque"], figure=1)
-plot(bouncingSphere, ["result.positionVector", "result.normalVector", "result.forceVector", "result.torqueVector"], figure=2)
+plot(bouncingSphere, ["result.penetration", "result.penetrationVelocity", "result.tangentialVelocity", "result.angularVelocity"], figure=1)
+plot(bouncingSphere, ["result.normalForce", "result.tangentialForce", "result.torque"], figure=2)
+plot(bouncingSphere, ["result.positionVector", "result.normalVector", "result.forceVector", "result.torqueVector"], figure=3)
 
 end