Add new feature to fluid2d: EMS

core/fluid2d.py

@@ -296,6 +296,7 @@ def signal_handler(signal, frame):
             # check for blow-up
             if self.model.diags['maxspeed'] > 1e3:
                 self.stop = True
+                self.blow_up = True
                 if self.myrank == 0:
                     print()
                     print('max|u| > 1000, blow-up detected, stopping')

core/output.py

@@ -26,7 +26,6 @@ def __init__(self, param, grid, diag, flxlist=None):
         self.template = self.expdir+'/%s_his_%03i.nc'
         if self.nbproc > 1:
             self.hisfile = self.template % (self.expname, self.myrank)
-            self.hisfile_joined = '%s/%s_his.nc' % (self.expdir, self.expname)
         else:
             self.hisfile = '%s/%s_his.nc' % (self.expdir, self.expname)
 
@@ -43,7 +42,6 @@ def __init__(self, param, grid, diag, flxlist=None):
             template = self.expdir+'/%s_flx_%03i.nc'
             if self.nbproc > 1:
                 self.flxfile = template % (self.expname, self.myrank)
-                self.flxfile_joined = '%s/%s_flx.nc' % (self.expdir, self.expname)
             else:
                 self.flxfile = '%s/%s_flx.nc' % (self.expdir, self.expname)
 
@@ -154,13 +152,13 @@ def dump_diag(self):
     def join(self):
         if self.nbproc > 1:
             filename = self.hisfile.split('his')[0]+'his'
-            join(filename)
+            self.hisfile = join(filename)
             if self.diag_fluxes:
                 filename = self.flxfile.split('flx')[0]+'flx'
-                join(filename)
-                
-            
-            
+                self.flxfile = join(filename)
+
+
+
 class NcfileIO(object):
     """Allow to create() and write() a Netcdf file of 'history' type,
     i.e. a set of model 2D snapshots. Variables were originally the
@@ -186,7 +184,11 @@ def create(self):
             # store all the parameters as NetCDF global attributes
             dparam = self.param.__dict__
             for k in dparam.keys():
-                if type(dparam[k]) == type([0, 1]):
+                if k == "ems":
+                    # The EMS stores itself and does not need to be stored in the netCDF file
+                    pass
+
+                elif type(dparam[k]) == type([0, 1]):
                     # it is not straightforward to store a list in a netcdf file
                     pass
 
@@ -250,7 +252,8 @@ def join(filename):
     ''' Join history files without having to mpirun
 
     Useful when the run has been broken and one wants to join
-    things from an interactive session '''
+    things from an interactive session.
+    Return the name of the new, joined history file. '''
 
     template = filename+'_%03i.nc'
     hisfile_joined = filename+'.nc'
@@ -390,3 +393,5 @@ def join(filename):
         os.remove(ncfile)
 
     print('-'*50)
+
+    return hisfile_joined

core/param.py

@@ -3,6 +3,9 @@
 import sys
 import getopt
 
+# Local import
+import ems
+
 
 class Param(object):
     """class to set up the default parameters value of the model
@@ -17,10 +20,16 @@ class Param(object):
 
     """
 
-    def __init__(self, defaultfile):
-        """defaultfile is a sequel, it's no longer used the default file is
-        systematically the defaults.json located in the fluid2d/core
+    def __init__(self, defaultfile=None, ems_file=""):
+        """Load default parameters and optionally experiment parameters.
+
+        The parameter `defaultfile` is no longer used and exists only
+        for backwards compatibility.  The default file is always
+        the file `defaults.json` located in the core-folder of fluid2d.
 
+        The parameter `ems_file` takes optionally the name of an
+        experiment file.  If given, the Experiment Management System
+        (EMS) is activated and the experiment file is parsed.
         """
 
         import grid
@@ -42,6 +51,11 @@ def __init__(self, defaultfile):
         else:
             self.print_param = False
 
+        if ems_file:
+            self.ems = ems.EMS(ems_file)
+        else:
+            self.ems = None
+
     def set_parameters(self, namelist):
         avail = {}
         doc = {}
@@ -76,6 +90,11 @@ def manall(self):
             self.man(p)
 
     def checkall(self):
+        if self.ems:
+            # Only create a new database entry once, not by every core
+            if self.myrank == 0:
+                self.ems.initialize(self.datadir)
+            self.expname = self.ems.get_expname()
         for p, avail in self.avail.items():
             if getattr(self, p) in avail:
                 # the parameter 'p' is well set
@@ -109,6 +128,36 @@ def copy(self, obj, list_param):
                 missing.append(k)
         return missing
 
+    def get_experiment_parameters(self):
+        """Return the experiment parameters dictionary loaded by the EMS.
+
+        The EMS must be activated in the constructor of `Param` to use
+        this method.  It returns the dictionary of experiment parameters
+        and exits.  It is advised to use, when possible, the method
+        `loop_experiment_parameters` instead.
+        """
+        return self.ems.parameters
+
+    def loop_experiment_parameters(self):
+        """Iterate over the experiment parameters loaded by the EMS.
+
+        In every iteration, this method returns a new dictionary of
+        experiment parameters containing a combination of the values
+        specified in the experiment file.  This experiment file for the
+        EMS must be specified in the constructor of `Param`.  If only
+        one value is given for every parameter in the experiment file,
+        the method `get_experiment_parameters` can be used instead.
+        The ID of the experiment is increased in every iteration.
+        """
+        while self.ems.parameters:
+            yield self.ems.parameters
+            self.ems.setup_next_parameters()
+
+    def finalize(self, fluid2d):
+        """Invoke the finalize method of the EMS if activated."""
+        if self.ems:
+            self.ems.finalize(fluid2d)
+
 
 if __name__ == "__main__":
     param = Param('default.xml')

core/plotting.py

@@ -32,6 +32,15 @@ def __init__(self, param, grid, var, diag):
                 # Define a dummy cax to prevent create_fig from failing
                 self.cax = [None, None]
 
+        if param.npx * param.npy > 1 and self.generate_mp4:
+            if self.myrank == 0:
+                print(
+                    'Warning: It is not possible to generate an mp4-file when '
+                    'fluid2d runs on multiple cores.\n'
+                    'The parameter generate_mp4 is automatically changed to False.'
+                )
+            self.generate_mp4 = False
+
         nh = self.nh
         nx = param.nx
         ny = param.ny
@@ -219,9 +228,10 @@ def update_fig(self, t, dt, kt):
             self.process.stdin.write(string)
 
     def finalize(self):
-        """ do nothing but close the mp4 thread if any"""
+        """Close the mp4 thread if any and the current figure."""
         if self.generate_mp4:
             self.process.communicate()
+        plt.close()
 
     def set_cax(self, z):
         """ set self.cax, the color range"""

core/restart.py

@@ -78,16 +78,10 @@ def __init__(self, param, grid, f2d, launch=True):
         f2d.output.template = self.expdir +'/%s_%02i_his' % (
             self.expname, self.nextrestart)+'_%03i.nc'
         f2d.output.hisfile = f2d.output.template % (self.myrank)
-        f2d.output.hisfile_joined = self.expdir + '/%s_%02i_his.nc' % (
-            self.expname, self.nextrestart)
         if self.diag_fluxes:
             f2d.output.template = self.expdir +'/%s_%02i_flx' % (
                 self.expname, self.nextrestart)+'_%03i.nc'
             f2d.output.flxfile = f2d.output.template % (self.myrank)
-            f2d.output.flxfile_joined = self.expdir + '/%s_%02i_flx.nc' % (
-                self.expname, self.nextrestart)
-
-
 
         # split the integration in 'ninterrestart' intervals and
         # save a restart at the end of each

ems/EMShellExtensions.py

@@ -0,0 +1,42 @@
+"""Extensions for the EMShell of fluid2d
+
+The EMShell is the command-line interface of the fluid2d Experiment
+Management System (EMS).  Its functionality for data analysis can be
+extended by adding new functions to the dictionary "extra_tools" in the
+EMShell.  These extensions are defined here, but they can also be
+imported from other files.
+
+Author: Markus REINERT, June 2019
+"""
+
+import numpy as np
+import netCDF4 as nc
+from scipy.fftpack import fft, fftshift, fftfreq
+
+
+def get_strongest_wavenumber_y(his_filename: str) -> float:
+    """Calculate the most intense wavenumber in y-direction.
+
+    This function opens the given history-file, performs a Fourier
+    transform in y on the masked streamfunction psi and returns the
+    highest wavenumber which has at some point in time the highest
+    intensity apart from the wavenumber zero."""
+    # Open history file and load the data
+    dataset_his = nc.Dataset(his_filename)
+    ny = dataset_his.ny
+    dy = dataset_his.Ly / ny
+    # Save the data as a masked numpy array
+    psi = dataset_his["psi"][:]
+    psi.mask = 1 - dataset_his["msk"][:]
+    # Set length of zero-padded signal (use ny for no zero-padding)
+    fft_ny = ny
+    # Caculate zero-padded Fourier-transform in y
+    fft_psi = fftshift(fft(psi, n=fft_ny, axis=1), axes=1)
+    # Its sampling frequency is dky = 1/dy/fft_ny
+    # Calculate the corresponding axis in Fourier-space
+    ky = fftshift(fftfreq(fft_ny, dy))
+    # Remove the zero-frequency because it is always very large
+    fft_psi[:, ky==0] = 0
+    # Calculate the frequency of maximal intensity (apart from the zero frequency)
+    ky_max = np.abs(ky)[np.argmax(np.max(np.abs(fft_psi), axis=2), axis=1)]
+    return np.max(ky_max)