Merge lp:~kahosato93/fluidity/firedrake-mesh into lp:fluidity/firedrake

[Kaho Sato]

Mesh files created along the creation of a mesh instance are now cached.

[lp:~kahosato93/fluidity/firedrake-mesh]

4328. By Kaho Sato <email address hidden> on 2013-09-10

Now mesh class uses cached mesh file if it exists

[Florian Rathgeber (florian-rathgeber)]

Excellent! No objections, merging...

review: Approve

=== modified file 'python/firedrake/mesh.py'
--- python/firedrake/mesh.py	2013-09-09 16:01:01 +0000
+++ python/firedrake/mesh.py	2013-09-10 15:49:10 +0000
@@ -7,57 +7,80 @@
 import firedrake
 
 
-class _get_msh_file(object):
-
-    """Class that can be used as a context manager of .geo and .msh files,
-    that removes these files as soon as they are no longer useful.
-    Given a source code and dimension of the mesh,
-    returns the object representing the temporary .msh file.
+_cachedir = os.path.join(tempfile.gettempdir(),
+                         'firedrake-mesh-cache-uid%d' % os.getuid())
+if not os.path.exists(_cachedir):
+    os.makedirs(_cachedir)
+
+_exts = [".node", ".ele", ".face"]
+_pexts = _exts + [".halo"]
+
+
+def _get_msh_file(source, name, dimension):
+
+    """Given a source code, name and dimension  of the mesh,
+    returns the name of the file that contains necessary information to build
+    a mesh class. The mesh class would call _from_file method on this file
+    to contruct itself.
     """
 
-    def __init__(self, source, dimension):
-        self._source = source
-        self._dimension = dimension
-
-    def __enter__(self):
-        if MPI.comm.rank == 0:
-            self._input = tempfile.NamedTemporaryFile('w+r', suffix='.geo')
-            self._input.write(self._source)
-            self._input.flush()
-            self._input.seek(0)
-            self._output = tempfile.NamedTemporaryFile('w+r', suffix='.msh')
-            result = self._output.name
-            subprocess.call(['gmsh', self._input.name,
-                            "-"+str(self._dimension), '-o', self._output.name])
+    if MPI.comm.rank == 0:
+        input = os.path.join(_cachedir, name+'.geo')
+        if not os.path.exists(input):
+            with open(input, 'w') as f:
+                f.write(source)
+
+        output = os.path.join(_cachedir, name)
+
+        if not os.path.exists(output+'.msh'):
+
+            #Writing of the output file.
             if MPI.parallel:
                 from mpi4py import MPI as _MPI
-                _MPI.COMM_SELF.Spawn('gmsh', args=[self._input.name,
-                        "-"+str(self._dimension), '-o', self._output.name])
+                _MPI.COMM_SELF.Spawn('gmsh', args=[input,
+                                     "-" + str(dimension), '-o', output+'.msh'])
+            else:
+                subprocess.call(['gmsh', input, "-"+str(dimension),
+                                 '-o', output+'.msh'])
+
+        if MPI.parallel:
+            if not _triangled(output, _exts):
                 gmsh2triangle = os.path.split(firedrake.__file__)[0] +\
                     "/../../bin/gmsh2triangle"
-                subprocess.call([gmsh2triangle, self._output.name])
-                result = os.path.splitext(self._output.name)[0][5:]
+                subprocess.call([gmsh2triangle, output+'.msh'], cwd=_cachedir)
+
+            basename = output + "_" + str(MPI.comm.size)
+            #Deal with decomposition.
+            #fldecomp would always name the decomposed triangle files
+            #in a same way.(meshname_rank.node, rather than
+            #meshname_size_rank.node).
+            #To go around this without creating triangle files everytime,
+            #we can make a simlink meshname_size.node which points to
+            #the file meshname.node.
+            for ext in _exts:
+                if not os.path.lexists(basename + ext):
+                    os.symlink(output+ext, basename+ext)
+            if not all([_triangled(basename + '_' + str(r), _pexts)
+                        for r in xrange(MPI.comm.size)]):
                 fldecomp = os.path.split(firedrake.__file__)[0] +\
                     "/../../bin/fldecomp"
                 subprocess.call([fldecomp, '-n', str(MPI.comm.size), '-m',
-                        'triangle', result])
-                result += ".node"
-                MPI.comm.bcast(result, root=0)
-
-        else:
-            result = MPI.comm.bcast(None, root=0)
-        return result
-
-    def __exit__(self, type, value, traceback):
-        if MPI.comm.rank == 0:
-            if MPI.parallel:
-                import glob
-                files = glob.glob(os.path.splitext(self._output.name)[0][5:] +
-                        '*')
-                for file in files:
-                    os.remove(file)
-            self._input.close()
-            self._output.close()
+                                 'triangle', basename])
+
+            output = basename + ".node"
+            MPI.comm.bcast(output, root=0)
+
+    #Not processor-0
+    else:
+        output = MPI.comm.bcast(None, root=0)
+
+    return output if MPI.parallel else output + '.msh'
+
+
+def _triangled(basename, exts):
+    """ Checks if the mesh of the given basename has already been decomposed.
+    """
+    return all(map(lambda ext: os.path.exists(basename + ext), exts))
 
 
 class UnitSquareMesh(Mesh):
@@ -82,9 +105,10 @@
                 Line{1}; Layers{%d};
             }
             """ % (nx, ny)
+        name = "unitsquare_%d_%d" % (nx, ny)
 
-        with _get_msh_file(source, 2) as output:
-            super(UnitSquareMesh, self).__init__(output)
+        output = _get_msh_file(source, name, 2)
+        super(UnitSquareMesh, self).__init__(output)
 
 
 class UnitCubeMesh(Mesh):
@@ -112,9 +136,10 @@
                 Surface{5}; Layers{%d};
             }
                        """ % (nx, ny, nz)
+        name = "unitcube_%d_%d_%d" % (nx, ny, nz)
 
-        with _get_msh_file(source, 3) as output:
-            super(UnitCubeMesh, self).__init__(output)
+        output = _get_msh_file(source, name, 3)
+        super(UnitCubeMesh, self).__init__(output)
 
 
 class UnitCircleMesh(Mesh):
@@ -135,9 +160,10 @@
                     Line{1};Layers{%d};
             }
             """ % (2*resolution, resolution * 4)
+        name = "unitcircle_%d" % resolution
 
-        with _get_msh_file(source, 2) as output:
-            super(UnitCircleMesh, self).__init__(output)
+        output = _get_msh_file(source, name, 2)
+        super(UnitCircleMesh, self).__init__(output)
 
 
 class UnitIntervalMesh(Mesh):
@@ -178,8 +204,10 @@
             Surface Loop(16) = {12, 14, 10, 8};
             Volume(16) = {16};
             """
-        with _get_msh_file(source, 3) as output:
-            super(UnitTetrahedronMesh, self).__init__(output)
+        name = "unittetra"
+
+        output = _get_msh_file(source, name, 3)
+        super(UnitTetrahedronMesh, self).__init__(output)
 
 
 class UnitTriangleMesh(Mesh):
@@ -196,5 +224,6 @@
             Line Loop(5) = {1, 2, 3};
             Plane Surface(5) = {5};
             """
-        with _get_msh_file(source, 2) as output:
-            super(UnitTriangleMesh, self).__init__(output)
+        name = "unittri"
+        output = _get_msh_file(source, name, 2)
+        super(UnitTriangleMesh, self).__init__(output)