From 54e4943353829708cf09d97ac1e3d085a75d8c81 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 30 Sep 2020 07:10:15 +0200 Subject: [PATCH 1/4] get rid of shell scripts --- processing/pre/geom_fromMinimalSurface.py | 21 ++------------ python/damask/_geom.py | 35 +++++++++++++++++++++++ python/tests/test_Geom.py | 18 ++++++++++++ 3 files changed, 55 insertions(+), 19 deletions(-) diff --git a/processing/pre/geom_fromMinimalSurface.py b/processing/pre/geom_fromMinimalSurface.py index 83d1a3684..bf0e2b754 100755 --- a/processing/pre/geom_fromMinimalSurface.py +++ b/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], - options.periods*2.0*np.pi*(np.arange(options.grid[1])+0.5)/options.grid[1], - 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:])]) +geom=damask.Geom.from_minimal_surface(options.grid,options.size,options.type,options.threshold, + options.periods,options.microstructure) damask.util.croak(geom) geom.save_ASCII(sys.stdout if name is None else name,compress=False) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 6a7bc3c0e..aaa534df0 100644 --- a/python/damask/_geom.py +++ b/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. diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 21eb38cdd..27829ac34 100644 --- a/python/tests/test_Geom.py +++ b/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) From a59e64a8e459b779a439edb0fc4022084d1de1d5 Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Wed, 30 Sep 2020 17:27:49 -0400 Subject: [PATCH 2/4] renamed TPMS and added more from additional references --- python/damask/_geom.py | 154 +++++++++++++++++++++++++------------- python/tests/test_Geom.py | 36 +++++++-- 2 files changed, 134 insertions(+), 56 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index aaa534df0..e760e40c3 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -120,6 +120,29 @@ class Geom: return np.unique(self.material).size + @staticmethod + def load(fname): + """ + Read a VTK rectilinear grid. + + Parameters + ---------- + fname : str or or pathlib.Path + Geometry file to read. + Valid extension is .vtr, it will be appended if not given. + + """ + v = VTK.load(fname if str(fname).endswith('.vtr') else str(fname)+'.vtr') + comments = v.get_comments() + grid = np.array(v.vtk_data.GetDimensions())-1 + bbox = np.array(v.vtk_data.GetBounds()).reshape(3,2).T + + return Geom(material = v.get('material').reshape(grid,order='F'), + size = bbox[1] - bbox[0], + origin = bbox[0], + comments=comments) + + @staticmethod def load_ASCII(fname): """ @@ -184,29 +207,6 @@ class Geom: return Geom(material.reshape(grid,order='F'),size,origin,comments) - @staticmethod - def load(fname): - """ - Read a VTK rectilinear grid. - - Parameters - ---------- - fname : str or or pathlib.Path - Geometry file to read. - Valid extension is .vtr, it will be appended if not given. - - """ - v = VTK.load(fname if str(fname).endswith('.vtr') else str(fname)+'.vtr') - comments = v.get_comments() - grid = np.array(v.vtk_data.GetDimensions())-1 - bbox = np.array(v.vtk_data.GetBounds()).reshape(3,2).T - - return Geom(material = v.get('material').reshape(grid,order='F'), - size = bbox[1] - bbox[0], - origin = bbox[0], - comments=comments) - - @staticmethod def _find_closest_seed(seeds, weights, point): return np.argmin(np.sum((np.broadcast_to(point,(len(seeds),3))-seeds)**2,axis=1) - weights) @@ -292,10 +292,52 @@ class Geom: ) + _minimal_surface = \ + {'Schwarz P': lambda x,y,z: np.cos(x) + np.cos(y) + np.cos(z), + 'Double Primitive': lambda x,y,z: ( 0.5 * (np.cos(x)*np.cos(y) + np.cos(y)*np.cos(z) + np.cos(z)*np.cos(x)) + + 0.2 * (np.cos(2*x) + np.cos(2*y) + np.cos(2*z)) ), + 'Schwarz D': lambda x,y,z: ( np.sin(x)*np.sin(y)*np.sin(z) + + np.sin(x)*np.cos(y)*np.cos(z) + + np.cos(x)*np.cos(y)*np.sin(z) + + np.cos(x)*np.sin(y)*np.cos(z) ), + 'Complementary D': lambda x,y,z: ( np.cos(3*x+y)*np.cos(z) - np.sin(3*x-y)*np.sin(z) + np.cos(x+3*y)*np.cos(z) + + np.sin(x-3*y)*np.sin(z) + np.cos(x-y)*np.cos(3*z) - np.sin(x+y)*np.sin(3*z) ), + 'Double Diamond': lambda x,y,z: 0.5 * (np.sin(x)*np.sin(y) + + np.sin(y)*np.sin(z) + + np.sin(z)*np.sin(x) + + np.cos(x) * np.cos(y) * np.cos(z) ), + 'Dprime': lambda x,y,z: 0.5 * ( np.cos(x)*np.cos(y)*np.cos(z) + + np.cos(x)*np.sin(y)*np.sin(z) + + np.sin(x)*np.cos(y)*np.sin(z) + + np.sin(x)*np.sin(y)*np.cos(z) + - np.sin(2*x)*np.sin(2*y) + - np.sin(2*y)*np.sin(2*z) + - np.sin(2*z)*np.sin(2*x) ) - 0.2, + 'Gyroid': lambda x,y,z: np.cos(x)*np.sin(y) + np.cos(y)*np.sin(z) + np.cos(z)*np.sin(x), + 'Gprime': lambda x,y,z : ( np.sin(2*x)*np.cos(y)*np.sin(z) + + np.sin(2*y)*np.cos(z)*np.sin(x) + + np.sin(2*z)*np.cos(x)*np.sin(y) ) + 0.32, + 'Karcher K': lambda x,y,z: ( 0.3 * ( np.cos(x) + np.cos(y) + np.cos(z) + + np.cos(x)*np.cos(y) + np.cos(y)*np.cos(z) + np.cos(z)*np.cos(x) ) + - 0.4 * ( np.cos(2*x) + np.cos(2*y) + np.cos(2*z) ) ) + 0.2, + 'Lidinoid': lambda x,y,z: 0.5 * ( np.sin(2*x)*np.cos(y)*np.sin(z) + + np.sin(2*y)*np.cos(z)*np.sin(x) + + np.sin(2*z)*np.cos(x)*np.sin(y) + - np.cos(2*x)*np.cos(2*y) + - np.cos(2*y)*np.cos(2*z) + - np.cos(2*z)*np.cos(2*x) ) + 0.15, + 'Neovius': lambda x,y,z: ( 3 * (np.cos(x)+np.cos(y)+np.cos(z)) + + 4 * np.cos(x)*np.cos(y)*np.cos(z) ), + 'Fisher-Koch S': lambda x,y,z: ( np.cos(2*x)*np.sin( y)*np.cos( z) + + np.cos( x)*np.cos(2*y)*np.sin( z) + + np.sin( x)*np.cos( y)*np.cos(2*z) ), + } + + @staticmethod def from_minimal_surface(grid,size,surface,threshold=0.0,periods=1,materials=(1,2)): """ - Generate geometry from definition of minimal surface. + Generate geometry from definition of triply periodic minimal surface. Parameters ---------- @@ -303,7 +345,7 @@ class Geom: 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'} + surface : {'primitive', 'gyroid', 'lidinoid', 'neovius', 'diamond', 'doublediamond'} Type of the minimal surface. threshold : float, optional. Threshold of the minimal surface. Defaults to 0.0. @@ -312,21 +354,53 @@ class Geom: materials : (int, int), optional Material IDs. Defaults to (1,2). + References + ---------- + Surface curvature in triply-periodic minimal surface architectures as + a distinct design parameter in preparing advanced tissue engineering scaffolds + Sébastien B G Blanquer, Maike Werner, Markus Hannula, Shahriar Sharifi, + Guillaume P R Lajoinie, David Eglin, Jari Hyttinen, André A Poot, and Dirk W Grijpma + 10.1088/1758-5090/aa6553 + + Triply Periodic Bicontinuous Cubic Microdomain Morphologies by Symmetries + Meinhard Wohlgemuth, Nataliya Yufa, James Hoffman, and Edwin L. Thomas + 10.1021/ma0019499 + + Minisurf – A minimal surface generator for finite element modeling and additive manufacturing + Meng-Ting Hsieh, Lorenzo Valdevit + 10.1016/j.simpa.2020.100026 + """ - 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]), + return Geom(material = np.where(threshold < Geom._minimal_surface[surface](x,y,z),materials[1],materials[0]), size = size, comments = util.execution_stamp('Geom','from_minimal_surface'), ) + def save(self,fname,compress=True): + """ + Generates vtk rectilinear grid. + + Parameters + ---------- + fname : str, optional + Filename to write. If no file is given, a string is returned. + Valid extension is .vtr, it will be appended if not given. + compress : bool, optional + Compress with zlib algorithm. Defaults to True. + + """ + v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) + v.add(self.material.flatten(order='F'),'material') + v.add_comments(self.comments) + + v.save(fname if str(fname).endswith('.vtr') else str(fname)+'.vtr',parallel=False,compress=compress) + + def save_ASCII(self,fname,compress=None): """ Writes a geom file. @@ -399,26 +473,6 @@ class Geom: f.write(f'{reps} of {former}\n') - def save(self,fname,compress=True): - """ - Generates vtk rectilinear grid. - - Parameters - ---------- - fname : str, optional - Filename to write. If no file is given, a string is returned. - Valid extension is .vtr, it will be appended if not given. - compress : bool, optional - Compress with zlib algorithm. Defaults to True. - - """ - v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) - v.add(self.material.flatten(order='F'),'material') - v.add_comments(self.comments) - - v.save(fname if str(fname).endswith('.vtr') else str(fname)+'.vtr',parallel=False,compress=compress) - - def show(self): """Show on screen.""" v = VTK.from_rectilinearGrid(self.grid,self.size,self.origin) diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index 27829ac34..884ae72e7 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -362,7 +362,19 @@ class TestGeom: assert np.all(geom.material == material) - @pytest.mark.parametrize('surface',['primitive','gyroid','diamond']) + @pytest.mark.parametrize('surface',['Schwarz P', + 'Double Primitive', + 'Schwarz D', + 'Complementary D', + 'Double Diamond', + 'Dprime', + 'Gyroid', + 'Gprime', + 'Karcher K', + 'Lidinoid', + 'Neovius', + 'Fisher-Koch S', + ]) def test_minimal_surface_basic_properties(self,surface): grid = np.random.randint(60,100,3) size = np.ones(3)+np.random.rand(3) @@ -373,8 +385,20 @@ class TestGeom: 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) + @pytest.mark.parametrize('surface,threshold',[('Schwarz P',0), + ('Double Primitive',-1./6.), + ('Schwarz D',0), + ('Complementary D',0), + ('Double Diamond',-0.133), + ('Dprime',-0.0395), + ('Gyroid',0), + ('Gprime',0.22913), + ('Karcher K',0.17045), + ('Lidinoid',0.14455), + ('Neovius',0), + ('Fisher-Koch S',0), + ]) + def test_minimal_surface_volume(self,surface,threshold): + grid = np.ones(3,dtype=int)*64 + geom = Geom.from_minimal_surface(grid,np.ones(3),surface,threshold) + assert np.isclose(np.count_nonzero(geom.material==1)/np.prod(geom.grid),.5,rtol=1e-3) From f2048087954aaaf33f15f946e67653342e7d4c0d Mon Sep 17 00:00:00 2001 From: Philip Eisenlohr Date: Wed, 30 Sep 2020 17:45:02 -0400 Subject: [PATCH 3/4] adapted geom_fromMinimalSurface.py to new TPMS names --- processing/pre/geom_fromMinimalSurface.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/processing/pre/geom_fromMinimalSurface.py b/processing/pre/geom_fromMinimalSurface.py index bf0e2b754..2b40940b7 100755 --- a/processing/pre/geom_fromMinimalSurface.py +++ b/processing/pre/geom_fromMinimalSurface.py @@ -11,7 +11,7 @@ scriptName = os.path.splitext(os.path.basename(__file__))[0] scriptID = ' '.join([scriptName,damask.version]) -minimal_surfaces = ['primitive','gyroid','diamond'] +minimal_surfaces = list(damask.Geom._minimal_surface.keys()) # -------------------------------------------------------------------- # MAIN From b29f22f51398b7c36608ab5b693239e00c19cb41 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Thu, 1 Oct 2020 09:25:32 +0200 Subject: [PATCH 4/4] documenting the actually available TMPS --- python/damask/_geom.py | 20 ++++++++++++++++++-- 1 file changed, 18 insertions(+), 2 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index e760e40c3..2ccbb1988 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -345,14 +345,30 @@ class Geom: 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', 'lidinoid', 'neovius', 'diamond', 'doublediamond'} - Type of the minimal surface. + surface : str + Type of the minimal surface. See notes for details. 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). + + Notes + ----- + The following triply-periodic minimal surfaces are implemented: + - Schwarz P + - Double Primitive + - Schwarz D + - Complementary D + - Double Diamond + - Dprime + - Gyroid + - Gprime + - Karcher K + - Lidinoid + - Neovius + - Fisher-Koch S References ----------