angle computation updated

lcharleux · Jan 11, 2018 · bc872a0 · bc872a0
1 parent 876d90f
commit bc872a0
Show file tree

Hide file tree

Showing 3 changed files with 27 additions and 20 deletions.
diff --git a/argiope/mesh.py b/argiope/mesh.py
@@ -62,6 +62,11 @@ def get_angles(self):
                            for i in range(len(self.faces))])
   angles = property(get_angles)  
 
+  def get_optimal_angles(self):
+    return  np.concatenate([ELEMENTS[e].optimal_angles 
+            for e in self.faces_types])
+  optimal_angles = property(get_optimal_angles)  
+
   def get_surfaces(self):
     return self.faces
   surfaces = property(get_surfaces)
@@ -97,7 +102,7 @@ def get_surfaces(self):
                                 "tri3",
                                 "tri3",
                                 "tri3"]),
-        simplices = np.array([[0, 1, 3, 4]])),
+        simplices = np.array([[0, 1, 2, 3]])),
     "pyra5": Element3D(
         nvert = 5,
         edges = np.array([[0, 1], [1, 2], [2, 3], [3, 0],
@@ -383,7 +388,7 @@ def centroids_and_volumes(self, sort_index = True):
       elif etype_info.space == 3:          
         simplices_volumes =  (np.cross(edges[:,:,0], 
                                        edges[:,:,1], axis = 2) 
-                             * edges[:,:, 2]).sum(axis = 2)
+                             * edges[:,:, 2]).sum(axis = 2) / 6.
       elements_volumes = simplices_volumes.sum(axis = 1)
       elements_centroids = ((simplices_volumes.reshape(*simplices_volumes.shape, 1) 
                           * simplices_centroids).sum(axis = 1) 
@@ -401,7 +406,7 @@ def centroids_and_volumes(self, sort_index = True):
     if sort_index: out.sort_index(inplace = True)
     return out.sort_index(axis= 1)
 
-  def angles(self):
+  def angles(self, zfill = 3):
     """
     Returns the internal angles of all elements and the associated statistics 
     """
@@ -430,11 +435,13 @@ def angles(self):
       angles_df = pd.DataFrame(index = index, 
                                data = angles, 
                                columns = pd.MultiIndex.from_product(
-      [["angles"], ["a{0}".format(s) for s in range(angles_info.shape[0])]]))
+      [["angles"], ["a" + "{0}".format(s).zfill(zfill) 
+              for s in range(angles_info.shape[0])]]))
       deviation_df = pd.DataFrame(index = index, 
                                data = deviation, 
                                columns = pd.MultiIndex.from_product(
-      [["deviation"], ["d{0}".format(s) for s in range(angles_info.shape[0])]]))
+      [["deviation"], ["d" + "{0}".format(s).zfill(zfill) 
+              for s in range(angles_info.shape[0])]]))
 
       df = pd.concat([angles_df, deviation_df], axis = 1).sort_index(axis = 1)
       df["stats", "max_angle"] = df.angles.max(axis = 1)

diff --git a/tests/mesh_example.py b/tests/mesh_example.py
@@ -7,17 +7,17 @@
                     [2., 0., 0.], #3
                     [1., 1., 0.], #4
                     [0., 1., 0.], #5
-                    [0., 0., -1.], #6
-                    [1., 0., -1.], #7 
-                    [2., 0., -1.], #8
-                    [1., 1., -1.], #9
-                    [0., 1., -1.], #10
-                    [0., 0., -2.], #11
-                    [1., 0., -2.], #12
-                    [2., 0., -2.], #13
-                    [1., 1., -2.], #14
-                    [0., 1., -2.], #15
-                    [1., 0., -3.], #16
+                    [0., 0., 1.], #6
+                    [1., 0., 1.], #7 
+                    [2., 0., 1.], #8
+                    [1., 1., 1.], #9
+                    [0., 1., 1.], #10
+                    [0., 0., 2.], #11
+                    [1., 0., 2.], #12
+                    [2., 0., 2.], #13
+                    [1., 1., 2.], #14
+                    [0., 1., 2.], #15
+                    [1., 0., 3.], #16
                     ])
 
 # NODE LABELS
@@ -32,7 +32,7 @@
          [2, 3, 4, 0, 0, 0, 0, 0], #2 = TRI3
          [6, 7, 9, 10, 11, 12, 14, 15], # 3 = HEXA8
          [7, 8, 9, 12, 13, 14, 0, 0], # 4 = PRISM6
-         [12, 13, 16, 14, 0, 0, 0, 0], # 5 = TETRA4
+         [12, 13, 14, 16, 0, 0, 0, 0], # 5 = TETRA4
         ]
 
 elabels = np.arange(1, len(conn) + 1)

diff --git a/tests/test_mesh.py b/tests/test_mesh.py
@@ -45,9 +45,9 @@ def test_volumes(self):
     Tests if volume computations are correct.
     """     
     vout = mesh.centroids_and_volumes().volume.values
-    vin = np.array([1., .5])
-    self.assertEqual( 
-       ((vin -vout)**2).sum(), 0 ) 
+    vin = np.array([1., .5, 1., .5, 1./6.])
+    self.assertTrue( 
+       ((vin -vout)**2).sum() < 1.e10, 0 ) 
 
 if __name__ == '__main__':
     unittest.main()