minor

Author: liminchen

friction/FrictionEnergy.py

@@ -3,23 +3,23 @@
 
 epsv = 1e-3
 
-def f0(vhatnorm, epsv, hhat):
-    if vhatnorm >= epsv:
-        return vhatnorm * hhat
+def f0(vbarnorm, epsv, hhat):
+    if vbarnorm >= epsv:
+        return vbarnorm * hhat
     else:
-        vhatnormhhat = vhatnorm * hhat
+        vbarnormhhat = vbarnorm * hhat
         epsvhhat = epsv * hhat
-        return vhatnormhhat * vhatnormhhat * (-vhatnormhhat / 3.0 + epsvhhat) / (epsvhhat * epsvhhat) + epsvhhat / 3.0
+        return vbarnormhhat * vbarnormhhat * (-vbarnormhhat / 3.0 + epsvhhat) / (epsvhhat * epsvhhat) + epsvhhat / 3.0
 
-def f1_div_vhatnorm(vhatnorm, epsv):
-    if vhatnorm >= epsv:
-        return 1.0 / vhatnorm
+def f1_div_vbarnorm(vbarnorm, epsv):
+    if vbarnorm >= epsv:
+        return 1.0 / vbarnorm
     else:
-        return (-vhatnorm + 2.0 * epsv) / (epsv * epsv)
+        return (-vbarnorm + 2.0 * epsv) / (epsv * epsv)
 
-def f_hess_term(vhatnorm, epsv):
-    if vhatnorm >= epsv:
-        return -1.0 / (vhatnorm * vhatnorm)
+def f_hess_term(vbarnorm, epsv):
+    if vbarnorm >= epsv:
+        return -1.0 / (vbarnorm * vbarnorm)
     else:
         return -1.0 / (epsv * epsv)
 
@@ -28,29 +28,29 @@ def val(v, mu_lambda, hhat, n):
     T = np.identity(2) - np.outer(n, n) # tangent of slope is constant
     for i in range(0, len(v)):
         if mu_lambda[i] > 0:
-            vhat = np.transpose(T).dot(v[i])
-            sum += mu_lambda[i] * f0(np.linalg.norm(vhat), epsv, hhat)
+            vbar = np.transpose(T).dot(v[i])
+            sum += mu_lambda[i] * f0(np.linalg.norm(vbar), epsv, hhat)
     return sum
 
 def grad(v, mu_lambda, hhat, n):
     g = np.array([[0.0, 0.0]] * len(v))
     T = np.identity(2) - np.outer(n, n) # tangent of slope is constant
     for i in range(0, len(v)):
         if mu_lambda[i] > 0:
-            vhat = np.transpose(T).dot(v[i])
-            g[i] = mu_lambda[i] * f1_div_vhatnorm(np.linalg.norm(vhat), epsv) * T.dot(vhat)
+            vbar = np.transpose(T).dot(v[i])
+            g[i] = mu_lambda[i] * f1_div_vbarnorm(np.linalg.norm(vbar), epsv) * T.dot(vbar)
     return g
 
 def hess(v, mu_lambda, hhat, n):
     IJV = [[0] * 0, [0] * 0, np.array([0.0] * 0)]
     T = np.identity(2) - np.outer(n, n) # tangent of slope is constant
     for i in range(0, len(v)):
         if mu_lambda[i] > 0:
-            vhat = np.transpose(T).dot(v[i])
-            vhatnorm = np.linalg.norm(vhat)
-            inner_term = f1_div_vhatnorm(vhatnorm, epsv) * np.identity(2)
-            if vhatnorm != 0:
-                inner_term += f_hess_term(vhatnorm, epsv) / vhatnorm * np.outer(vhat, vhat)
+            vbar = np.transpose(T).dot(v[i])
+            vbarnorm = np.linalg.norm(vbar)
+            inner_term = f1_div_vbarnorm(vbarnorm, epsv) * np.identity(2)
+            if vbarnorm != 0:
+                inner_term += f_hess_term(vbarnorm, epsv) / vbarnorm * np.outer(vbar, vbar)
             local_hess = mu_lambda[i] * T.dot(utils.make_PD(inner_term)).dot(np.transpose(T)) / hhat
             for c in range(0, 2):
                 for r in range(0, 2):

friction/simulator.py

@@ -15,6 +15,7 @@
 h = 0.01        # time step size in s
 DBC = []        # no nodes need to be fixed
 ground_n = np.array([0.1, 1.0])     # normal of the slope
+ground_n /= np.linalg.norm(ground_n)    # normalize ground normal vector just in case
 ground_o = np.array([0.0, -1.0])    # a point on the slope  
 mu = 0.11        # friction coefficient of the slope
 
@@ -33,7 +34,6 @@
 for i in DBC:
     is_DBC[i] = True
 contact_area = [side_len / n_seg] * len(x)     # perimeter split to each node
-ground_n /= np.linalg.norm(ground_n)    # normalize ground normal vector just in case
 
 # simulation with visualization
 resolution = np.array([900, 900])
@@ -56,7 +56,7 @@ def screen_projection(x):
     # fill the background and draw the square
     screen.fill((255, 255, 255))
     pygame.draw.aaline(screen, (0, 0, 255), screen_projection([ground_o[0] - 3.0 * ground_n[1], ground_o[1] + 3.0 * ground_n[0]]), 
-        screen_projection([ground_o[0] + 3.0 * ground_n[1], ground_o[1] - 3.0 * ground_n[0]]))   # ground
+        screen_projection([ground_o[0] + 3.0 * ground_n[1], ground_o[1] - 3.0 * ground_n[0]]))   # slope
     for eI in e:
         pygame.draw.aaline(screen, (0, 0, 255), screen_projection(x[eI[0]]), screen_projection(x[eI[1]]))
     for xI in x: