From a0e0f28e51db223110cb87d0ebdb99c95804b341 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sat, 8 Aug 2020 19:42:34 +0200 Subject: [PATCH] migrating shell scripts to library --- processing/pre/geom_vicinityOffset.py | 32 ++------- python/damask/_geom.py | 97 ++++++++++++++++++++------- python/tests/test_Geom.py | 41 ++++++++++- 3 files changed, 117 insertions(+), 53 deletions(-) diff --git a/processing/pre/geom_vicinityOffset.py b/processing/pre/geom_vicinityOffset.py index e30779d31..650e0093a 100755 --- a/processing/pre/geom_vicinityOffset.py +++ b/processing/pre/geom_vicinityOffset.py @@ -5,7 +5,6 @@ import sys from io import StringIO from optparse import OptionParser -from scipy import ndimage import numpy as np import damask @@ -15,18 +14,6 @@ scriptName = os.path.splitext(os.path.basename(__file__))[0] scriptID = ' '.join([scriptName,damask.version]) -def taintedNeighborhood(stencil,trigger=[],size=1): - - me = stencil[stencil.shape[0]//2] - if len(trigger) == 0: - return np.any(stencil != me) - if me in trigger: - trigger = set(trigger) - trigger.remove(me) - trigger = list(trigger) - return np.any(np.in1d(stencil,np.array(trigger))) - - #-------------------------------------------------------------------------------------------------- # MAIN #-------------------------------------------------------------------------------------------------- @@ -67,19 +54,12 @@ options.trigger = np.array(options.trigger, dtype=int) if filenames == []: filenames = [None] for name in filenames: - damask.util.report(scriptName,name) + damask.util.report(scriptName,name) - geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name) + geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name) - offset = np.nanmax(geom.microstructure) if options.offset is None else options.offset + damask.util.croak(geom.vicinity_offset(options.vicinity,options.offset,options.trigger, + True if options.mode is 'wrap' else False)) + geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:])) - damask.util.croak(geom.update(np.where(ndimage.filters.generic_filter( - geom.microstructure, - taintedNeighborhood, - size=1+2*options.vicinity,mode=options.mode, - extra_arguments=(), - extra_keywords={"trigger":options.trigger,"size":1+2*options.vicinity}), - geom.microstructure + offset,geom.microstructure))) - geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:])) - - geom.to_file(sys.stdout if name is None else name,pack=False) + geom.to_file(sys.stdout if name is None else name,pack=False) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 144a80b70..4416be401 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -278,7 +278,7 @@ class Geom: Parameters ---------- fname : str or file handle - geometry file to read. + Geometry file to read. """ try: @@ -345,15 +345,15 @@ class Geom: Parameters ---------- grid : numpy.ndarray of shape (3) - number of grid points in x,y,z direction. + Number of grid points in x,y,z direction. size : list or numpy.ndarray of shape (3) - physical size of the microstructure in meter. + Physical size of the microstructure in meter. seeds : numpy.ndarray of shape (:,3) - position of the seed points in meter. All points need to lay within the box. + Position of the seed points in meter. All points need to lay within the box. weights : numpy.ndarray of shape (seeds.shape[0]) - weights of the seeds. Setting all weights to 1.0 gives a standard Voronoi tessellation. + Weights of the seeds. Setting all weights to 1.0 gives a standard Voronoi tessellation. periodic : Boolean, optional - perform a periodic tessellation. Defaults to True. + Perform a periodic tessellation. Defaults to True. """ if periodic: @@ -391,13 +391,13 @@ class Geom: Parameters ---------- grid : numpy.ndarray of shape (3) - number of grid points in x,y,z direction. + Number of grid points in x,y,z direction. size : list or numpy.ndarray of shape (3) - physical size of the microstructure in meter. + Physical size of the microstructure in meter. seeds : numpy.ndarray of shape (:,3) - position of the seed points in meter. All points need to lay within the box. + Position of the seed points in meter. All points need to lay within the box. periodic : Boolean, optional - perform a periodic tessellation. Defaults to True. + Perform a periodic tessellation. Defaults to True. """ coords = grid_filters.cell_coord0(grid,size).reshape(-1,3) @@ -415,9 +415,9 @@ class Geom: Parameters ---------- fname : str or file handle - geometry file to write. + Geometry file to write. pack : bool, optional - compress geometry with 'x of y' and 'a to b'. + Compress geometry with 'x of y' and 'a to b'. """ header = self.get_header() @@ -481,7 +481,7 @@ class Geom: Parameters ---------- fname : str, optional - vtk file to write. If no file is given, a string is returned. + Vtk file to write. If no file is given, a string is returned. """ v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) @@ -508,9 +508,9 @@ class Geom: Parameters ---------- directions : iterable containing str - direction(s) along which the microstructure is mirrored. Valid entries are 'x', 'y', 'z'. + Direction(s) along which the microstructure is mirrored. Valid entries are 'x', 'y', 'z'. reflect : bool, optional - reflect (include) outermost layers. + Reflect (include) outermost layers. """ valid = {'x','y','z'} @@ -540,7 +540,7 @@ class Geom: Parameters ---------- grid : numpy.ndarray of shape (3) - number of grid points in x,y,z direction. + Number of grid points in x,y,z direction. """ #ToDo: self.add_comments('geom.py:scale v{}'.format(version) @@ -563,7 +563,7 @@ class Geom: Parameters ---------- stencil : int, optional - size of smoothing stencil. + Size of smoothing stencil. """ def mostFrequent(arr): @@ -596,9 +596,9 @@ class Geom: Parameters ---------- R : damask.Rotation - rotation to apply to the microstructure. + Rotation to apply to the microstructure. fill : int or float, optional - microstructure index to fill the corners. Defaults to microstructure.max() + 1. + Microstructure index to fill the corners. Defaults to microstructure.max() + 1. """ if fill is None: fill = np.nanmax(self.microstructure) + 1 @@ -631,11 +631,11 @@ class Geom: Parameters ---------- grid : numpy.ndarray of shape (3) - number of grid points in x,y,z direction. + Number of grid points in x,y,z direction. offset : numpy.ndarray of shape (3) - offset (measured in grid points) from old to new microstructue[0,0,0]. + Offset (measured in grid points) from old to new microstructure[0,0,0]. fill : int or float, optional - microstructure index to fill the corners. Defaults to microstructure.max() + 1. + Microstructure index to fill the corners. Defaults to microstructure.max() + 1. """ if fill is None: fill = np.nanmax(self.microstructure) + 1 @@ -658,14 +658,14 @@ class Geom: def substitute(self,from_microstructure,to_microstructure): """ - Substitude microstructure indices. + Substitute microstructure indices. Parameters ---------- from_microstructure : iterable of ints - microstructure indices to be substituted. + Microstructure indices to be substituted. to_microstructure : iterable of ints - new microstructure indices. + New microstructure indices. """ substituted = self.get_microstructure() @@ -674,3 +674,50 @@ class Geom: #ToDo: self.add_comments('geom.py:substitute v{}'.format(version) return self.update(substituted) + + + def vicinity_offset(self,vicinity=1,offset=None,trigger=[],periodic=True): + """ + Offset microstructure index of points in the vicinity of xxx. + + Different from themselves (or listed as triggers) within a given (cubic) vicinity, + i.e. within the region close to a grain/phase boundary. + ToDo: use include/exclude as in seeds.from_geom + + Parameters + ---------- + vicinity : int, optional + Voxel distance checked for presence of other microstructure. + Defaults to 1. + offset : int, optional + Offset (positive or negative) to tag microstructure indices, + defaults to microstructure.max() + 1. + trigger : list of ints, optional + List of microstructure indices triggering a change. + Defaults to [], meaning that different neigboors trigger a change. + periodic : Boolean, optional + Assume geometry to be periodic. Defaults to True. + + """ + def tainted_neighborhood(stencil,trigger): + + me = stencil[stencil.shape[0]//2] + if len(trigger) == 0: + return np.any(stencil != me) + if me in trigger: + trigger = set(trigger) + trigger.remove(me) + trigger = list(trigger) + return np.any(np.in1d(stencil,np.array(trigger))) + + offset_ = np.nanmax(self.microstructure) if offset is None else offset + mask = ndimage.filters.generic_filter(self.microstructure, + tainted_neighborhood, + size=1+2*vicinity, + mode=('wrap' if periodic else 'nearest'), + extra_keywords={'trigger':trigger}) + microstructure = np.ma.MaskedArray(self.microstructure + offset_, np.logical_not(mask)) + + #ToDo: self.add_comments('geom.py:vicinity_offset v{}'.format(version) + return self.update(microstructure) + diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index f5ffbd6bf..1cb29e0b0 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -77,9 +77,19 @@ class TestGeom: with pytest.raises(ValueError): default.update(default.microstructure[1:,1:,1:],size=np.ones(2)) - def test_invalid_microstructure(self,default): + def test_invalid_origin(self,default): with pytest.raises(ValueError): - default.update(default.microstructure[1]) + default.update(default.microstructure[1:,1:,1:],origin=np.ones(4)) + + def test_invalid_microstructure_size(self,default): + microstructure=np.ones((3,3)) + with pytest.raises(ValueError): + default.update(microstructure) + + def test_invalid_microstructure_type(self,default): + microstructure=np.random.randint(1,300,(3,4,5))==1 + with pytest.raises(TypeError): + default.update(microstructure) def test_invalid_homogenization(self,default): with pytest.raises(TypeError): @@ -171,7 +181,34 @@ class TestGeom: modified.canvas(modified.grid + grid_add) e = default.grid assert np.all(modified.microstructure[:e[0],:e[1],:e[2]] == default.microstructure) + + @pytest.mark.parametrize('trigger',[[1],[]]) + def test_vicinity_offset(self,trigger): + offset = np.random.randint(2,4) + vicinity = np.random.randint(2,4) + + g=np.random.randint(28,40,(3)) + m=np.ones(g,'i') + x=(g*np.random.permutation(np.array([.5,1,1]))).astype('i') + m[slice(0,x[0]),slice(0,x[1]),slice(0,x[2])]=2 + m2=copy.deepcopy(m) + for i in [0,1,2]: + m2[(np.roll(m,+vicinity,i)-m)!=0] +=offset + m2[(np.roll(m,-vicinity,i)-m)!=0] +=offset + if len(trigger) > 0: + m2[m==1]=1 + + geom = Geom(m,np.random.rand(3)) + geom.vicinity_offset(vicinity,offset,trigger=trigger) + + assert np.all(m2==geom.microstructure) + @pytest.mark.parametrize('periodic',[True,False]) + def test_vicinity_offset_invariant(self,default,periodic): + old = default.get_microstructure() + default.vicinity_offset(trigger=[old.max()+1,old.min()-1]) + assert np.all(old==default.microstructure) + @pytest.mark.parametrize('periodic',[True,False]) def test_tessellation_approaches(self,periodic): grid = np.random.randint(10,20,3)