improved 2dm file format compat.

pvelissariou1 · Jan 25, 2020 · 82d4dcc · 82d4dcc
1 parent b168219
commit 82d4dcc
Show file tree

Hide file tree

Showing 6 changed files with 69 additions and 41 deletions.
diff --git a/pyschism/mesh/base.py b/pyschism/mesh/base.py
@@ -317,14 +317,14 @@ def sms2dm(self):
     def tria3(self):
         return np.array(
             [list(map(self.get_node_index, index))
-             for id, index in self._triangles.items()])
+             for index in self._triangles.values()])
 
     @property
     @lru_cache
     def quad4(self):
         return np.array(
             [list(map(self.get_node_index, index))
-             for id, index in self._quads.items()])
+             for index in self._quads.values()])
 
     @property
     def logger(self):
@@ -419,9 +419,9 @@ def _sms2dm(self):
         if self.crs is not None:
             description += f"; {self.crs.srs}"
         return {
-            "nodes": self.nodes,
-            "elements": self.elements,
-            "description": description,
+            "ND": self._nodes,
+            "E3T": self._triangles,
+            "E4Q": self._quads,
         }
 
     @_coords.setter

diff --git a/pyschism/mesh/gmesh.py b/pyschism/mesh/gmesh.py
@@ -15,7 +15,9 @@
 
 
 class Gmesh(EuclideanMesh2D):
-
+    """
+    boundaries = {ibtype: {id: {'indexes': [i0, ..., in], 'properties': object }}
+    """
     def __init__(
         self,
         coords,
@@ -25,7 +27,7 @@ def __init__(
         crs=None,
         description=None,
         boundaries=None,
-    ): 
+    ):
         super().__init__(coords, triangles, quads, values, crs, description)
         self._boundaries = boundaries
 
@@ -62,7 +64,6 @@ def delete_boundary_data(self, ibtype, id):
 
     def write_boundaries(self, path, overwrite=False):
         path = pathlib.Path(path)
-        # print(path)
         if path.exists() and not overwrite:
             msg = "Destination path exists and overwrite=False"
             raise IOError(msg)
@@ -203,7 +204,9 @@ def _grd(self):
     @property
     def _sms2dm(self):
         sms2dm = super()._sms2dm
-        sms2dm.update({'boundaries': self.boundaries})
+        def nodestrings(boundaries):
+            return
+        sms2dm.update({'nodestrings': nodestrings(self.boundaries)})
         return sms2dm
 
     @property
@@ -212,9 +215,6 @@ def _boundaries(self):
 
     @_boundaries.setter
     def _boundaries(self, boundaries):
-        """
-        boundaries = {ibtype: {id: {'indexes': [i0, ..., in], 'properties': object }}
-        """
         self.clear_boundaries() # clear
         if boundaries is not None:
             for ibtype, bnds in boundaries.items():

diff --git a/pyschism/mesh/sms2dm.py b/pyschism/mesh/sms2dm.py
@@ -35,6 +35,25 @@ def string(sms2dm):
     f += boundaries(sms2dm)
     return f
 
+
+def graph(sms2dm):
+    f = "MESH2D\n"
+    # TODO: Make faster using np.array2string
+    f += triangular_elements(sms2dm)
+    f += quadrilateral_elements(sms2dm)
+    f += nodes(sms2dm)
+    return f
+
+def nodes(sms2dm):
+    assert all(int(id) > 0 for id in sms2dm['ND'])
+    f = ''
+    for id, (coords, value)in sms2dm['ND'].items():
+        f += f"ND {int(id)} {id:d}"
+        f += f"{coords[0]:<.16E} "
+        f += f"{coords[1]:<.16E} "
+        f += f"{-value:<.16E}\n"
+    return f
+
 def boundaries(sms2dm):
     f = ''
     if 'boundaries' in sms2dm.keys():
@@ -43,13 +62,14 @@ def boundaries(sms2dm):
                 f += nodestring(bnd['indexes'])
     return f
 
-def geom_string(geom_type, geom):
+def geom_string(geom_type, sms2dm):
     assert geom_type in ['E3T', 'E4Q', 'E6T', 'E8Q', 'E9Q']
+    assert all(int(id) > 0 for id in sms2dm[geom_type])
     f = ''
-    for i in range(len(geom)):
-        f += f"{geom_type} {i + 1} "
-        for j in range(len(geom[i, :])):
-            f += f"{geom[i, j]+1} "
+    for id, geom in sms2dm[geom_type].items():
+        f += f"{geom_type} {id} "
+        for j in range(len(geom)):
+            f += f"{geom[j]} "
         f += "\n"
     return f
 
@@ -67,23 +87,5 @@ def triangular_elements(geom):
 def quadrilateral_elements(geom):
     f = ''
     if geom is not None:
-        for i in range(len(geom)):
-            f += f"E4Q {i + 1} "
-            for j in range(len(geom[i, :])):
-                f += f"{geom[i, j]+1} "
-            f += "\n"
+        f += geom_string("E4Q", geom)
     return f
-
-def graph(sms2dm):
-    f = "MESH2D\n"
-    # TODO: Make faster using np.array2string
-    if 'triangles' in sms2dm:
-        f += geom_string("E3T", sms2dm['triangles'])
-    if 'quads' in sms2dm:
-        f += geom_string("E4Q", sms2dm['quads'])
-    for i in range(len(sms2dm['coords'])):
-        f += f"ND {i + 1} "
-        f += f"{sms2dm['coords'][i][0]:<.16E} "
-        f += f"{sms2dm['coords'][i][1]:<.16E} "
-        f += f"{-sms2dm['values'][i]:<.16E}\n"
-    return f
diff --git a/tests/coverage.svg b/tests/coverage.svg
diff --git a/tests/mesh/test_base.py b/tests/mesh/test_base.py
@@ -205,8 +205,10 @@ def test_open_fmt_grd(self):
         msh = EuclideanMesh2D.open(gr3.absolute(), fmt='grd')
         self.assertIsInstance(msh, EuclideanMesh2D)
 
-    def test_open_fmt_2dm(self):
-        pass
+    def test_sms2dm(self):
+        m = EuclideanMesh2D(self.coords, self.triangles, self.quads)
+        self.assertIsInstance(m.sms2dm, str)
+
 
 
 if __name__ == '__main__':

diff --git a/tests/mesh/test_gmesh.py b/tests/mesh/test_gmesh.py
@@ -287,7 +287,7 @@ def test_write_boundaries_raises(self):
         msh.logger.debug('coverage')
         self.assertRaises(IOError, msh.write_boundaries, shp)
 
-    def test_gr3(self):
+    def test_grd(self):
         boundaries = {
             None: {
                 0: {
@@ -310,5 +310,29 @@ def test_gr3(self):
             boundaries=boundaries)
         self.assertIsInstance(msh.grd, str)
 
+    def test_sms2dm(self):
+        boundaries = {
+            None: {
+                0: {
+                    'indexes':
+                        [139510, 140443, 140461, 140462, 141993, 150761]
+                    }
+                },
+            0: {
+                0: {
+                    'indexes':
+                        [139510, 140443, 140461, 140462, 141993, 150761],
+                    'properties': {'fake_velocites': None},
+                    }
+                }
+            }
+        msh = Gmesh(
+            self.coords,
+            self.triangles,
+            self.quads,
+            crs="EPSG:3395",
+            boundaries=boundaries)
+        self.assertIsInstance(msh.sms2dm, str)
+
 if __name__ == '__main__':
     unittest.main()