From e61c1a027b2d326c595baf5d0de5d41a0936e00b Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Sun, 29 Mar 2020 19:12:23 +0200 Subject: [PATCH] avoid detour via shell --- python/damask/_geom.py | 97 +++++++++++++++++++++++++++++++++++---- python/tests/test_Geom.py | 21 +++++++++ 2 files changed, 109 insertions(+), 9 deletions(-) diff --git a/python/damask/_geom.py b/python/damask/_geom.py index 6bc92e300..797b9dc81 100644 --- a/python/damask/_geom.py +++ b/python/damask/_geom.py @@ -1,11 +1,15 @@ import sys from io import StringIO +import multiprocessing +from functools import partial import numpy as np -from scipy import ndimage +from scipy import ndimage,spatial from . import VTK from . import util +from . import Environment +from . import grid_filters class Geom: @@ -50,7 +54,7 @@ class Geom: def update(self,microstructure=None,size=None,origin=None,rescale=False): """ - Updates microstructure and size. + Update microstructure and size. Parameters ---------- @@ -113,7 +117,7 @@ class Geom: def set_comments(self,comments): """ - Replaces all existing comments. + Replace all existing comments. Parameters ---------- @@ -127,7 +131,7 @@ class Geom: def add_comments(self,comments): """ - Appends comments to existing comments. + Append comments to existing comments. Parameters ---------- @@ -140,7 +144,7 @@ class Geom: def set_microstructure(self,microstructure): """ - Replaces the existing microstructure representation. + Replace the existing microstructure representation. Parameters ---------- @@ -159,7 +163,7 @@ class Geom: def set_size(self,size): """ - Replaces the existing size information. + Replace the existing size information. Parameters ---------- @@ -179,7 +183,7 @@ class Geom: def set_origin(self,origin): """ - Replaces the existing origin information. + Replace the existing origin information. Parameters ---------- @@ -196,7 +200,7 @@ class Geom: def set_homogenization(self,homogenization): """ - Replaces the existing homogenization index. + Replace the existing homogenization index. Parameters ---------- @@ -264,7 +268,7 @@ class Geom: @staticmethod def from_file(fname): """ - Reads a geom file. + Read a geom file. Parameters ---------- @@ -325,6 +329,81 @@ class Geom: return Geom(microstructure.reshape(grid),size,origin,homogenization,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) + @staticmethod + def from_Laguerre_tessellation(grid,size,seeds,weights,periodic=True): + """ + Generate geometry from Laguerre tessellation. + + Parameters + ---------- + grid : 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 microstructure in meter. + seeds : numpy.ndarray of shape (:,3) + 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. + periodic : Boolean, optional + perform a periodic tessellation. Defaults to True. + + """ + if periodic: + weights_p = np.tile(weights,27).flatten(order='F') # Laguerre weights (1,2,3,1,2,3,...,1,2,3) + seeds_p = np.vstack((seeds -np.array([size[0],0.,0.]),seeds, seeds +np.array([size[0],0.,0.]))) + seeds_p = np.vstack((seeds_p-np.array([0.,size[1],0.]),seeds_p,seeds_p+np.array([0.,size[1],0.]))) + seeds_p = np.vstack((seeds_p-np.array([0.,0.,size[2]]),seeds_p,seeds_p+np.array([0.,0.,size[2]]))) + coords = grid_filters.cell_coord0(grid*3,size*3,-size).reshape(-1,3,order='F') + + else: + weights_p = weights.flatten() + seeds_p = seeds + coords = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') + + pool = multiprocessing.Pool(processes = int(Environment().options['DAMASK_NUM_THREADS'])) + result = pool.map_async(partial(Geom._find_closest_seed,seeds_p,weights_p), [coord for coord in coords]) + pool.close() + pool.join() + microstructure = np.array(result.get()) + + if periodic: + microstructure = microstructure.reshape(grid[0]*3,grid[1]*3,grid[2]*3) + microstructure = microstructure[grid[0]:grid[0]*2,grid[1]:grid[1]*2,grid[2]:grid[2]*2]%seeds.shape[0] + else: + microstructure = microstructure.reshape(grid) + + #comments = 'geom.py:from_Laguerre_tessellation v{}'.format(version) + return Geom(microstructure+1,size,homogenization=1) + + + @staticmethod + def from_Voronoi_tessellation(grid,size,seeds,periodic=True): + """ + Generate geometry from Voronoi tessellation. + + Parameters + ---------- + grid : 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 microstructure in meter. + seeds : numpy.ndarray of shape (:,3) + 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. + + """ + coords = grid_filters.cell_coord0(grid,size).reshape(-1,3,order='F') + KDTree = spatial.cKDTree(seeds,boxsize=size) if periodic else spatial.cKDTree(seeds) + devNull,microstructure = KDTree.query(coords) + + #comments = 'geom.py:from_Voronoi_tessellation v{}'.format(version) + return Geom(microstructure.reshape(grid)+1,size,homogenization=1) + + def to_file(self,fname,pack=None): """ Writes a geom file. diff --git a/python/tests/test_Geom.py b/python/tests/test_Geom.py index a6fd57cb6..1b1529e80 100644 --- a/python/tests/test_Geom.py +++ b/python/tests/test_Geom.py @@ -97,3 +97,24 @@ class TestGeom: reference = os.path.join(reference_dir,'scale_{}.geom'.format(tag)) if update: modified.to_file(reference) assert geom_equal(modified,Geom.from_file(reference)) + + @pytest.mark.parametrize('periodic',[(True),(False)]) + def test_tessellation(self,periodic): + grid = np.random.randint(10,20,3) + size = np.random.random(3) + 1.0 + N_seeds= np.random.randint(10,30) + seeds = np.random.rand(N_seeds,3) * np.broadcast_to(size,(N_seeds,3)) + Voronoi = Geom.from_Voronoi_tessellation( grid,size,seeds, periodic) + Laguerre = Geom.from_Laguerre_tessellation(grid,size,seeds,np.ones(N_seeds),periodic) + assert geom_equal(Laguerre,Voronoi) + + def test_Laguerre(self): + grid = np.random.randint(10,20,3) + size = np.random.random(3) + 1.0 + N_seeds= np.random.randint(10,30) + seeds = np.random.rand(N_seeds,3) * np.broadcast_to(size,(N_seeds,3)) + weights= np.full((N_seeds),-np.inf) + ms = np.random.randint(1, N_seeds+1) + weights[ms-1] = np.random.random() + Laguerre = Geom.from_Laguerre_tessellation(grid,size,seeds,weights,np.random.random()>0.5) + assert np.all(Laguerre.microstructure == ms)