get rid of shell scripts

processing/pre/geom_fromMinimalSurface.py

@@ -4,8 +4,6 @@ import os
 import sys
 from optparse import OptionParser
 
-import numpy as np
-
 import damask
 
 
@@ -15,13 +13,6 @@ scriptID   = ' '.join([scriptName,damask.version])
 
 minimal_surfaces = ['primitive','gyroid','diamond']
 
-surface = {
-            'primitive': lambda x,y,z: np.cos(x)+np.cos(y)+np.cos(z),
-            'gyroid':    lambda x,y,z: np.sin(x)*np.cos(y)+np.sin(y)*np.cos(z)+np.cos(x)*np.sin(z),
-            'diamond':   lambda x,y,z: np.cos(x-y)*np.cos(z)+np.sin(x+y)*np.sin(z),
-          }
-
-
 # --------------------------------------------------------------------
 #                                MAIN
 # --------------------------------------------------------------------
@@ -71,16 +62,8 @@ parser.set_defaults(type = minimal_surfaces[0],
 name = None if filename == [] else filename[0]
 damask.util.report(scriptName,name)
 
-x,y,z = np.meshgrid(options.periods*2.0*np.pi*(np.arange(options.grid[0])+0.5)/options.grid[0],
+geom=damask.Geom.from_minimal_surface(options.grid,options.size,options.type,options.threshold,
-                    options.periods*2.0*np.pi*(np.arange(options.grid[1])+0.5)/options.grid[1],
+                                      options.periods,options.microstructure)
-                    options.periods*2.0*np.pi*(np.arange(options.grid[2])+0.5)/options.grid[2],
-                    indexing='xy',sparse=True)
-
-microstructure = np.where(options.threshold < surface[options.type](x,y,z),
-                          options.microstructure[1],options.microstructure[0])
-
-geom=damask.Geom(microstructure,options.size,
-                 comments=[scriptID + ' ' + ' '.join(sys.argv[1:])])
 damask.util.croak(geom)
 
 geom.save_ASCII(sys.stdout if name is None else name,compress=False)

python/damask/_geom.py

@@ -292,6 +292,41 @@ class Geom:
                    )
 
 
+    @staticmethod
+    def from_minimal_surface(grid,size,surface,threshold=0.0,periods=1,materials=(1,2)):
+        """
+        Generate geometry from definition of minimal surface.
+
+        Parameters
+        ----------
+        grid : int numpy.ndarray of shape (3)
+            Number of grid points in x,y,z direction.
+        size : list or numpy.ndarray of shape (3)
+            Physical size of the geometry in meter.
+        surface : {'primitive', 'gyroid', 'diamond'}
+            Type of the minimal surface.
+        threshold : float, optional.
+            Threshold of the minimal surface. Defaults to 0.0.
+        periods : integer, optional.
+            Number of periods per unit cell. Defaults to 1.
+        materials : (int, int), optional
+            Material IDs. Defaults to (1,2).
+
+        """
+        s = {'primitive': lambda x,y,z: np.cos(x)+np.cos(y)+np.cos(z),
+             'gyroid':    lambda x,y,z: np.sin(x)*np.cos(y)+np.sin(y)*np.cos(z)+np.cos(x)*np.sin(z),
+             'diamond':   lambda x,y,z: np.cos(x-y)*np.cos(z)+np.sin(x+y)*np.sin(z),
+          }
+        x,y,z = np.meshgrid(periods*2.0*np.pi*(np.arange(grid[0])+0.5)/grid[0],
+                            periods*2.0*np.pi*(np.arange(grid[1])+0.5)/grid[1],
+                            periods*2.0*np.pi*(np.arange(grid[2])+0.5)/grid[2],
+                            indexing='ij',sparse=True)
+        return Geom(material = np.where(threshold < s[surface](x,y,z),materials[1],materials[0]),
+                    size     = size,
+                    comments = util.execution_stamp('Geom','from_minimal_surface'),
+                   )
+
+
     def save_ASCII(self,fname,compress=None):
         """
         Writes a geom file.

python/tests/test_Geom.py

@@ -360,3 +360,21 @@ class TestGeom:
         elif approach == 'Voronoi':
             geom = Geom.from_Voronoi_tessellation(grid,size,seeds,            periodic=np.random.random()>0.5)
         assert np.all(geom.material == material)
+
+
+    @pytest.mark.parametrize('surface',['primitive','gyroid','diamond'])
+    def test_minimal_surface_basic_properties(self,surface):
+        grid = np.random.randint(60,100,3)
+        size = np.ones(3)+np.random.rand(3)
+        threshold = np.random.rand()-.5
+        periods = np.random.randint(2)+1
+        materials = np.random.randint(0,40,2)
+        geom = Geom.from_minimal_surface(grid,size,surface,threshold,periods,materials)
+        assert geom.material.max() == materials.max() and geom.material.min() == materials.min() \
+               and (geom.size == size).all() and (geom.grid == grid).all()
+
+    @pytest.mark.parametrize('surface',['primitive','gyroid','diamond'])
+    def test_minimal_surface_volume(self,surface):
+        grid = np.ones(3,dtype='i')*64
+        geom = Geom.from_minimal_surface(grid,np.ones(3),surface)
+        assert np.isclose(np.count_nonzero(geom.material==1)/np.prod(geom.grid),.5)