From 446ac03b072117ec8979973027df0667206084ed Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Sun, 23 Aug 2020 18:46:01 -0400 Subject: [PATCH] All geom methods are now out-of-place, i.e. return an updated duplicate (to allow for daisy chaining). * Added comments when methods acted. * Added diff method * Added flip method * Fixed add_primitive inversion bug (again...) * Fixed cell centering bug in add_primitive * Added missing tests --- python/damask/_geom.py | 255 ++++++++++++++++++++++++++++---------- python/tests/test_Geom.py | 124 ++++++++++-------- 2 files changed, 263 insertions(+), 116 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 0f1c35e42..c4bbbbdef 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -7,6 +7,7 @@ from functools import partial import numpy as np from scipy import ndimage,spatial +from . import version from . import environment from . import Rotation from . import VTK @@ -24,15 +25,15 @@ class Geom: Parameters ---------- microstructure : numpy.ndarray - microstructure array (3D) + Microstructure array (3D) size : list or numpy.ndarray - physical size of the microstructure in meter. + Physical size of the microstructure in meter. origin : list or numpy.ndarray, optional - physical origin of the microstructure in meter. + Physical origin of the microstructure in meter. homogenization : int, optional - homogenization index. + Homogenization index. comments : list of str, optional - comments lines. + Comment lines. """ self.set_microstructure(microstructure) @@ -50,7 +51,7 @@ class Geom: f'origin x y z: {util.srepr(self.get_origin()," ")}', f'# microstructures: {self.N_microstructure}', f'max microstructure: {np.nanmax(self.microstructure)}', - ]) + ]+self.get_comments()) def __copy__(self): @@ -63,20 +64,57 @@ class Geom: return self.__copy__() - def update(self,microstructure=None,size=None,origin=None,rescale=False): + def duplicate(self,microstructure=None,size=None,origin=None,comments=None,autosize=False): + """ + Create a duplicate having updated microstructure, size, and origin. + + Parameters + ---------- + microstructure : numpy.ndarray, optional + Microstructure array (3D). + size : list or numpy.ndarray, optional + Physical size of the microstructure in meter. + origin : list or numpy.ndarray, optional + Physical origin of the microstructure in meter. + comments : list of str, optional + Comment lines. + autosize : bool, optional + Ignore size parameter and rescale according to change of grid points. + + """ + if size is not None and autosize: + raise ValueError('Auto sizing conflicts with explicit size parameter.') + + grid_old = self.get_grid() + dup = self.copy() + dup.set_microstructure(microstructure) + dup.set_origin(origin) + + if comments is not None: + dup.set_comments(comments) + + if size is not None: + dup.set_size(size) + elif autosize: + dup.set_size(dup.get_grid()/grid_old*self.get_size()) + + return dup + + + def update(self,microstructure=None,size=None,origin=None,autosize=False): """ Update microstructure and size. Parameters ---------- microstructure : numpy.ndarray, optional - microstructure array (3D). + Microstructure array (3D). size : list or numpy.ndarray, optional - physical size of the microstructure in meter. + Physical size of the microstructure in meter. origin : list or numpy.ndarray, optional - physical origin of the microstructure in meter. - rescale : bool, optional - ignore size parameter and rescale according to change of grid points. + Physical origin of the microstructure in meter. + autosize : bool, optional + Ignore size parameter and rescale according to change of grid points. """ grid_old = self.get_grid() @@ -85,15 +123,15 @@ class Geom: unique_old = self.N_microstructure max_old = np.nanmax(self.microstructure) - if size is not None and rescale: - raise ValueError('Either set size explicitly or rescale automatically') + if size is not None and autosize: + raise ValueError('Auto sizing conflicts with explicit size parameter.') self.set_microstructure(microstructure) self.set_origin(origin) if size is not None: self.set_size(size) - elif rescale: + elif autosize: self.set_size(self.get_grid()/grid_old*self.size) message = [f'grid a b c: {util.srepr(grid_old," x ")}'] @@ -124,6 +162,40 @@ class Geom: return util.return_message(message) + def diff(self,other): + """ + Report property differences of self relative to other. + + Parameters + ---------- + other : Geom + Geometry to compare self against. + + """ + message = [] + if np.any(other.get_grid() != self.get_grid()): + message.append(util.delete(f'grid a b c: {util.srepr(other.get_grid()," x ")}')) + message.append(util.emph( f'grid a b c: {util.srepr( self.get_grid()," x ")}')) + + if np.any(other.get_size() != self.get_size()): + message.append(util.delete(f'size x y z: {util.srepr(other.get_size()," x ")}')) + message.append(util.emph( f'size x y z: {util.srepr( self.get_size()," x ")}')) + + if np.any(other.get_origin() != self.get_origin()): + message.append(util.delete(f'origin x y z: {util.srepr(other.get_origin()," ")}')) + message.append(util.emph( f'origin x y z: {util.srepr( self.get_origin()," ")}')) + + if other.N_microstructure != self.N_microstructure: + message.append(util.delete(f'# microstructures: {other.N_microstructure}')) + message.append(util.emph( f'# microstructures: { self.N_microstructure}')) + + if np.nanmax(other.microstructure) != np.nanmax(self.microstructure): + message.append(util.delete(f'max microstructure: {np.nanmax(other.microstructure)}')) + message.append(util.emph( f'max microstructure: {np.nanmax( self.microstructure)}')) + + return util.return_message(message) + + def set_comments(self,comments): """ Replace all existing comments. @@ -131,7 +203,7 @@ class Geom: Parameters ---------- comments : list of str - new comments. + All comments. """ self.comments = [] @@ -145,7 +217,7 @@ class Geom: Parameters ---------- comments : list of str - new comments. + New comments. """ self.comments += [str(c) for c in comments] if isinstance(comments,list) else [str(comments)] @@ -172,6 +244,10 @@ class Geom: else: self.microstructure = np.copy(microstructure) + if self.microstructure.dtype == float and \ + np.all(self.microstructure == self.microstructure.astype(int).astype(float)): + self.microstructure = self.microstructure.astype(int) + if len(self.microstructure.shape) != 3: raise ValueError(f'Invalid microstructure shape {microstructure.shape}') elif self.microstructure.dtype not in np.sctypes['float'] + np.sctypes['int']: @@ -185,7 +261,7 @@ class Geom: Parameters ---------- size : list or numpy.ndarray - physical size of the microstructure in meter. + Physical size of the microstructure in meter. """ if size is None: @@ -205,7 +281,7 @@ class Geom: Parameters ---------- origin : list or numpy.ndarray - physical origin of the microstructure in meter + Physical origin of the microstructure in meter. """ if origin is not None: @@ -222,12 +298,12 @@ class Geom: Parameters ---------- homogenization : int - homogenization index + Homogenization index. """ if homogenization is not None: if not isinstance(homogenization,int) or homogenization < 1: - raise TypeError(f'Invalid homogenization {homogenization}') + raise TypeError(f'Invalid homogenization {homogenization}.') else: self.homogenization = homogenization @@ -419,8 +495,11 @@ class Geom: else: microstructure = microstructure.reshape(grid) - #ToDo: comments = 'geom.py:from_Laguerre_tessellation v{}'.format(version) - return Geom(microstructure+1,size,homogenization=1) + return Geom(microstructure+1, + size, + homogenization=1, + comments=f'geom.py:from_Laguerre_tessellation v{version}', + ) @staticmethod @@ -444,8 +523,11 @@ class Geom: KDTree = spatial.cKDTree(seeds,boxsize=size) if periodic else spatial.cKDTree(seeds) devNull,microstructure = KDTree.query(coords) - #ToDo: comments = 'geom.py:from_Voronoi_tessellation v{}'.format(version) - return Geom(microstructure.reshape(grid)+1,size,homogenization=1) + return Geom(microstructure.reshape(grid)+1, + size, + homogenization=1, + comments=f'geom.py:from_Voronoi_tessellation v{version}', + ) def to_file(self,fname,pack=None): @@ -535,7 +617,7 @@ class Geom: def show(self): """Show raw content (as in file).""" - f=StringIO() + f = StringIO() self.to_file(f) f.seek(0) return ''.join(f.readlines()) @@ -565,10 +647,10 @@ class Geom: R : damask.Rotation, optional Rotation of primitive. Defaults to no rotation. inverse : Boolean, optional - Retain original microstructure within primitive and fill - outside. Defaults to False. + Retain original microstructure within primitive and fill outside. + Defaults to False. periodic : Boolean, optional - Repeat primitive over boundaries. Defaults to False. + Repeat primitive over boundaries. Defaults to True. """ # normalized 'radius' and center @@ -578,11 +660,11 @@ class Geom: (np.array(center) - self.origin)/self.size coords = grid_filters.cell_coord0(self.grid,np.ones(3)) \ - - (np.ones(3)*0.5 if periodic else c) # center if periodic + - ((np.ones(3)-(1./self.grid if np.array(center).dtype in np.sctypes['int'] else 0))*0.5 if periodic else c) # periodic center is always at CoG coords_rot = R.broadcast_to(tuple(self.grid))@coords with np.errstate(over='ignore',under='ignore'): - mask = np.where(np.sum(np.abs(coords_rot/r)**(2.0**exponent),axis=-1) < 1,True,False) + mask = np.where(np.sum(np.abs(coords_rot/r)**(2.0**exponent),axis=-1) <= 1.0,False,True) if periodic: # translate back to center mask = np.roll(mask,((c-np.ones(3)*.5)*self.grid).astype(int),(0,1,2)) @@ -590,7 +672,9 @@ class Geom: fill_ = np.full_like(self.microstructure,np.nanmax(self.microstructure)+1 if fill is None else fill) ms = np.ma.MaskedArray(fill_,np.logical_not(mask) if inverse else mask) - return self.update(ms) + return self.duplicate(ms, + comments=self.get_comments()+[f'geom.py:add_primitive v{version}'], + ) def mirror(self,directions,reflect=False): @@ -622,8 +706,41 @@ class Geom: if 'x' in directions: ms = np.concatenate([ms,ms[limits[0]:limits[1]:-1,:,:]],0) - #ToDo: self.add_comments('geom.py:mirror v{}'.format(version) - return self.update(ms,rescale=True) + return self.duplicate(ms, + comments=self.get_comments()+[f'geom.py:mirror {set(directions)} v{version}'], + autosize=True) + + + def flip(self,directions): + """ + Flip microstructure along given directions. + + Parameters + ---------- + directions : iterable containing str + Direction(s) along which the microstructure is flipped. + Valid entries are 'x', 'y', 'z'. + + """ + valid = {'x','y','z'} + if not all(isinstance(d, str) for d in directions): + raise TypeError('Directions are not of type str.') + elif not set(directions).issubset(valid): + raise ValueError(f'Invalid direction {set(directions).difference(valid)} specified.') + + limits = [None,None] + ms = self.get_microstructure() + + if 'z' in directions: + ms = ms[:,:,limits[0]:limits[1]:-1] + if 'y' in directions: + ms = ms[:,limits[0]:limits[1]:-1,:] + if 'x' in directions: + ms = ms[limits[0]:limits[1]:-1,:,:] + + return self.duplicate(ms, + comments=self.get_comments()+[f'geom.py:flip {set(directions)} v{version}'], + ) def scale(self,grid): @@ -636,17 +753,16 @@ class Geom: Number of grid points in x,y,z direction. """ - #ToDo: self.add_comments('geom.py:scale v{}'.format(version) - return self.update( - ndimage.interpolation.zoom( - self.microstructure, - grid/self.get_grid(), - output=self.microstructure.dtype, - order=0, - mode='nearest', - prefilter=False - ) - ) + return self.duplicate(ndimage.interpolation.zoom( + self.microstructure, + grid/self.get_grid(), + output=self.microstructure.dtype, + order=0, + mode='nearest', + prefilter=False + ), + comments=self.get_comments()+[f'geom.py:scale {grid} v{version}'], + ) def clean(self,stencil=3,mode='nearest',selection=None): @@ -672,15 +788,15 @@ class Geom: else: return me - #ToDo: self.add_comments('geom.py:clean v{}'.format(version) - return self.update(ndimage.filters.generic_filter( - self.microstructure, - mostFrequent, - size=(stencil if selection is None else stencil//2*2+1,)*3, - mode=mode, - extra_keywords=dict(selection=selection), - ).astype(self.microstructure.dtype) - ) + return self.duplicate(ndimage.filters.generic_filter( + self.microstructure, + mostFrequent, + size=(stencil if selection is None else stencil//2*2+1,)*3, + mode=mode, + extra_keywords=dict(selection=selection), + ).astype(self.microstructure.dtype), + comments=self.get_comments()+[f'geom.py:clean {stencil} {mode} v{version}'], + ) def renumber(self): @@ -689,8 +805,9 @@ class Geom: for i, oldID in enumerate(np.unique(self.microstructure)): renumbered = np.where(self.microstructure == oldID, i+1, renumbered) - #ToDo: self.add_comments('geom.py:renumber v{}'.format(version) - return self.update(renumbered) + return self.duplicate(renumbered, + comments=self.get_comments()+[f'geom.py:renumber v{version}'], + ) def rotate(self,R,fill=None): @@ -724,8 +841,11 @@ class Geom: origin = self.origin-(np.asarray(microstructure_in.shape)-self.grid)*.5 * self.size/self.grid - #ToDo: self.add_comments('geom.py:rotate v{}'.format(version) - return self.update(microstructure_in,origin=origin,rescale=True) + return self.duplicate(microstructure_in, + origin=origin, + comments=self.get_comments()+[f'geom.py:rotate {R.as_quaternion()} v{version}'], + autosize=True, + ) def canvas(self,grid=None,offset=None,fill=None): @@ -757,8 +877,11 @@ class Geom: canvas[ll[0]:ur[0],ll[1]:ur[1],ll[2]:ur[2]] = self.microstructure[LL[0]:UR[0],LL[1]:UR[1],LL[2]:UR[2]] - #ToDo: self.add_comments('geom.py:canvas v{}'.format(version) - return self.update(canvas,origin=self.origin+offset*self.size/self.grid,rescale=True) + return self.duplicate(canvas, + origin=self.origin+offset*self.size/self.grid, + comments=self.get_comments()+[f'geom.py:canvas {grid} {offset} {fill} v{version}'], + autosize=True, + ) def substitute(self,from_microstructure,to_microstructure): @@ -777,8 +900,9 @@ class Geom: for from_ms,to_ms in zip(from_microstructure,to_microstructure): substituted[self.microstructure==from_ms] = to_ms - #ToDo: self.add_comments('geom.py:substitute v{}'.format(version) - return self.update(substituted) + return self.duplicate(substituted, + comments=self.get_comments()+[f'geom.py:substitute v{version}'], + ) def vicinity_offset(self,vicinity=1,offset=None,trigger=[],periodic=True): @@ -819,9 +943,10 @@ class Geom: mask = ndimage.filters.generic_filter(self.microstructure, tainted_neighborhood, size=1+2*vicinity, - mode=('wrap' if periodic else 'nearest'), + 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) + return self.duplicate(microstructure, + comments=self.get_comments()+[f'geom.py:vicinity_offset {vicinity} v{version}'], + ) diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 4fb1ae00a..e568b4a53 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -28,6 +28,19 @@ def reference_dir(reference_dir_base): class TestGeom: + @pytest.mark.parametrize('flavor',['plain','explicit']) + def test_duplicate(self,default,flavor): + if flavor == 'plain': + modified = default.duplicate() + elif flavor == 'explicit': + modified = default.duplicate( + default.get_microstructure(), + default.get_size(), + default.get_origin() + ) + print(modified) + assert geom_equal(default,modified) + def test_update(self,default): modified = default.copy() modified.update( @@ -36,7 +49,7 @@ class TestGeom: default.get_origin() ) print(modified) - assert geom_equal(modified,default) + assert geom_equal(default,modified) @pytest.mark.parametrize('masked',[True,False]) def test_set_microstructure(self,default,masked): @@ -47,36 +60,38 @@ class TestGeom: def test_write_read_str(self,default,tmpdir): - default.to_file(str(tmpdir.join('default.geom'))) - new = Geom.from_file(str(tmpdir.join('default.geom'))) - assert geom_equal(new,default) + default.to_file(tmpdir/'default.geom') + new = Geom.from_file(tmpdir/'default.geom') + assert geom_equal(default,new) def test_write_read_file(self,default,tmpdir): - with open(tmpdir.join('default.geom'),'w') as f: + with open(tmpdir/'default.geom','w') as f: default.to_file(f) - with open(tmpdir.join('default.geom')) as f: + with open(tmpdir/'default.geom') as f: new = Geom.from_file(f) - assert geom_equal(new,default) + assert geom_equal(default,new) def test_write_show(self,default,tmpdir): - with open(tmpdir.join('str.geom'),'w') as f: + with open(tmpdir/'str.geom','w') as f: f.write(default.show()) - with open(tmpdir.join('str.geom')) as f: + with open(tmpdir/'str.geom') as f: new = Geom.from_file(f) - assert geom_equal(new,default) + assert geom_equal(default,new) + + def test_export_import_vtk(self,default,tmpdir): + default.to_vtk(tmpdir/'default') + assert geom_equal(default,Geom.from_vtk(tmpdir/'default.vtr')) - def test_export_vtk(self,default,tmpdir): - default.to_vtk(str(tmpdir.join('default'))) @pytest.mark.parametrize('pack',[True,False]) def test_pack(self,default,tmpdir,pack): - default.to_file(tmpdir.join('default.geom'),pack=pack) - new = Geom.from_file(tmpdir.join('default.geom')) + default.to_file(tmpdir/'default.geom',pack=pack) + new = Geom.from_file(tmpdir/'default.geom') assert geom_equal(new,default) def test_invalid_combination(self,default): with pytest.raises(ValueError): - default.update(default.microstructure[1:,1:,1:],size=np.ones(3), rescale=True) + default.update(default.microstructure[1:,1:,1:],size=np.ones(3), autosize=True) def test_invalid_size(self,default): with pytest.raises(ValueError): @@ -108,21 +123,36 @@ class TestGeom: ] ) def test_mirror(self,default,update,reference_dir,directions,reflect): - modified = default.copy() - modified.mirror(directions,reflect) + modified = default.mirror(directions,reflect) tag = f'directions={"-".join(directions)}_reflect={reflect}' reference = reference_dir/f'mirror_{tag}.geom' if update: modified.to_file(reference) - assert geom_equal(modified,Geom.from_file(reference)) + assert geom_equal(Geom.from_file(reference), + modified) + + @pytest.mark.parametrize('directions',[ + ['x'], + ['x','y','z'], + ['z','x','y'], + ['y','z'], + ] + ) + def test_flip(self,default,update,reference_dir,directions): + modified = default.flip(directions) + tag = f'directions={"-".join(directions)}' + reference = reference_dir/f'flip_{tag}.geom' + if update: modified.to_file(reference) + assert geom_equal(Geom.from_file(reference), + modified) @pytest.mark.parametrize('stencil',[1,2,3,4]) def test_clean(self,default,update,reference_dir,stencil): - modified = default.copy() - modified.clean(stencil) + modified = default.clean(stencil) tag = f'stencil={stencil}' reference = reference_dir/f'clean_{tag}.geom' if update: modified.to_file(reference) - assert geom_equal(modified,Geom.from_file(reference)) + assert geom_equal(Geom.from_file(reference), + modified) @pytest.mark.parametrize('grid',[ (10,11,10), @@ -134,33 +164,29 @@ class TestGeom: ] ) def test_scale(self,default,update,reference_dir,grid): - modified = default.copy() - modified.scale(grid) + modified = default.scale(grid) tag = f'grid={util.srepr(grid,"-")}' reference = reference_dir/f'scale_{tag}.geom' if update: modified.to_file(reference) - assert geom_equal(modified,Geom.from_file(reference)) + assert geom_equal(Geom.from_file(reference), + modified) def test_renumber(self,default): - modified = default.copy() - microstructure = modified.get_microstructure() + microstructure = default.get_microstructure() for m in np.unique(microstructure): microstructure[microstructure==m] = microstructure.max() + np.random.randint(1,30) - modified.update(microstructure) + modified = default.duplicate(microstructure) assert not geom_equal(modified,default) - modified.renumber() - assert geom_equal(modified,default) + assert geom_equal(default, + modified.renumber()) def test_substitute(self,default): - modified = default.copy() - microstructure = modified.get_microstructure() offset = np.random.randint(1,500) - microstructure += offset - modified.update(microstructure) + modified = default.duplicate(default.get_microstructure() + offset) assert not geom_equal(modified,default) - modified.substitute(np.arange(default.microstructure.max())+1+offset, - np.arange(default.microstructure.max())+1) - assert geom_equal(modified,default) + assert geom_equal(default, + modified.substitute(np.arange(default.microstructure.max())+1+offset, + np.arange(default.microstructure.max())+1)) @pytest.mark.parametrize('axis_angle',[np.array([1,0,0,86.7]), np.array([0,1,0,90.4]), np.array([0,0,1,90]), np.array([1,0,0,175]),np.array([0,-1,0,178]),np.array([0,0,1,180])]) @@ -168,38 +194,35 @@ class TestGeom: modified = default.copy() for i in range(np.rint(360/axis_angle[3]).astype(int)): modified.rotate(Rotation.from_axis_angle(axis_angle,degrees=True)) - assert geom_equal(modified,default) + assert geom_equal(default,modified) @pytest.mark.parametrize('Eulers',[[32.0,68.0,21.0], [0.0,32.0,240.0]]) def test_rotate(self,default,update,reference_dir,Eulers): - modified = default.copy() - modified.rotate(Rotation.from_Eulers(Eulers,degrees=True)) + modified = default.rotate(Rotation.from_Eulers(Eulers,degrees=True)) tag = f'Eulers={util.srepr(Eulers,"-")}' reference = reference_dir/f'rotate_{tag}.geom' if update: modified.to_file(reference) - assert geom_equal(modified,Geom.from_file(reference)) + assert geom_equal(Geom.from_file(reference), + modified) def test_canvas(self,default): + grid = default.grid grid_add = np.random.randint(0,30,(3)) - modified = default.copy() - modified.canvas(modified.grid + grid_add) - e = default.grid - assert np.all(modified.microstructure[:e[0],:e[1],:e[2]] == default.microstructure) + modified = default.canvas(grid + grid_add) + assert np.all(modified.microstructure[:grid[0],:grid[1],:grid[2]] == default.microstructure) @pytest.mark.parametrize('center1,center2',[(np.random.random(3)*.5,np.random.random(3)), (np.random.randint(4,8,(3)),np.random.randint(9,12,(3)))]) @pytest.mark.parametrize('diameter',[np.random.random(3)*.5, - np.random.randint(4,10,(3))]) + np.random.randint(4,10,(3))]) def test_add_primitive(self,diameter,center1,center2): """Same volume fraction for periodic microstructures and different center.""" o = np.random.random(3)-.5 g = np.random.randint(8,32,(3)) s = np.random.random(3)+.5 - G_1 = Geom(np.ones(g,'i'),s,o) - G_2 = Geom(np.ones(g,'i'),s,o) - G_1.add_primitive(diameter,center1,1) - G_2.add_primitive(diameter,center2,1) + G_1 = Geom(np.ones(g,'i'),s,o).add_primitive(diameter,center1,1) + G_2 = Geom(np.ones(g,'i'),s,o).add_primitive(diameter,center2,1) assert np.count_nonzero(G_1.microstructure!=2) == np.count_nonzero(G_2.microstructure!=2) @pytest.mark.parametrize('trigger',[[1],[]]) @@ -218,8 +241,7 @@ class TestGeom: if len(trigger) > 0: m2[m==1] = 1 - geom = Geom(m,np.random.rand(3)) - geom.vicinity_offset(vicinity,offset,trigger=trigger) + geom = Geom(m,np.random.rand(3)).vicinity_offset(vicinity,offset,trigger=trigger) assert np.all(m2==geom.microstructure)