DAMASK_EICMD/python/damask/seeds.py

155 lines
5.8 KiB
Python
Raw Normal View History

2020-11-16 05:42:23 +05:30
"""Functionality for generation of seed points for Voronoi or Laguerre tessellation."""
from typing import Tuple as _Tuple
2021-11-02 22:31:32 +05:30
from scipy import spatial as _spatial
2020-09-23 13:15:36 +05:30
import numpy as _np
from ._typehints import FloatSequence as _FloatSequence, IntSequence as _IntSequence, \
NumpyRngSeed as _NumpyRngSeed, IntCollection as _IntCollection
from . import util as _util
from . import grid_filters as _grid_filters
2020-09-23 13:15:36 +05:30
def from_random(size: _FloatSequence,
N_seeds: int,
cells: _IntSequence = None,
rng_seed: _NumpyRngSeed = None) -> _np.ndarray:
2020-09-23 13:15:36 +05:30
"""
Place seeds randomly in space.
2020-09-25 00:56:16 +05:30
2020-09-23 13:15:36 +05:30
Parameters
----------
size : sequence of float, len (3)
Physical size of the seeding domain.
2020-09-23 13:15:36 +05:30
N_seeds : int
Number of seeds.
cells : sequence of int, len (3), optional.
2020-12-08 05:06:41 +05:30
If given, ensures that each seed results in a grain when a standard Voronoi
tessellation is performed using the given grid resolution (i.e. size/cells).
rng_seed : {None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
2020-09-23 13:15:36 +05:30
A seed to initialize the BitGenerator. Defaults to None.
If None, then fresh, unpredictable entropy will be pulled from the OS.
2020-09-25 00:56:16 +05:30
Returns
-------
coords : numpy.ndarray, shape (N_seeds,3)
2021-04-24 21:30:57 +05:30
Seed coordinates in 3D space.
2020-09-23 13:15:36 +05:30
"""
size_ = _np.array(size,float)
rng = _np.random.default_rng(rng_seed)
2020-12-04 02:28:24 +05:30
if cells is None:
coords = rng.random((N_seeds,3)) * size_
2020-09-23 13:15:36 +05:30
else:
grid_coords = _grid_filters.coordinates0_point(cells,size).reshape(-1,3,order='F')
2020-12-04 02:28:24 +05:30
coords = grid_coords[rng.choice(_np.prod(cells),N_seeds, replace=False)] \
+ _np.broadcast_to(size_/_np.array(cells,int),(N_seeds,3))*(rng.random((N_seeds,3))*.5-.25) # wobble w/o leaving grid
2020-09-23 13:15:36 +05:30
2020-09-25 00:56:16 +05:30
return coords
2020-09-23 13:15:36 +05:30
def from_Poisson_disc(size: _FloatSequence,
N_seeds: int,
N_candidates: int,
distance: float,
periodic: bool = True,
rng_seed: _NumpyRngSeed = None) -> _np.ndarray:
2020-09-23 13:15:36 +05:30
"""
Place seeds according to a Poisson disc distribution.
2020-09-25 00:56:16 +05:30
2020-09-23 13:15:36 +05:30
Parameters
----------
size : sequence of float, len (3)
Physical size of the seeding domain.
2020-09-23 13:15:36 +05:30
N_seeds : int
Number of seeds.
N_candidates : int
Number of candidates to consider for finding best candidate.
distance : float
Minimum acceptable distance to other seeds.
2022-01-13 03:43:38 +05:30
periodic : bool, optional
Calculate minimum distance for periodically repeated grid.
rng_seed : {None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
2020-09-23 13:15:36 +05:30
A seed to initialize the BitGenerator. Defaults to None.
If None, then fresh, unpredictable entropy will be pulled from the OS.
2020-09-25 00:56:16 +05:30
Returns
-------
coords : numpy.ndarray, shape (N_seeds,3)
2021-04-24 21:30:57 +05:30
Seed coordinates in 3D space.
2020-09-23 13:15:36 +05:30
"""
rng = _np.random.default_rng(rng_seed)
2020-09-23 13:15:36 +05:30
coords = _np.empty((N_seeds,3))
coords[0] = rng.random(3) * _np.array(size,float)
2020-09-23 13:15:36 +05:30
2021-04-28 11:27:20 +05:30
s = 1
i = 0
progress = _util.ProgressBar(N_seeds+1,'',50)
2021-04-28 11:27:20 +05:30
while s < N_seeds:
2021-11-02 22:31:32 +05:30
i += 1
2020-09-23 13:15:36 +05:30
candidates = rng.random((N_candidates,3))*_np.broadcast_to(size,(N_candidates,3))
2021-04-28 11:27:20 +05:30
tree = _spatial.cKDTree(coords[:s],boxsize=size) if periodic else \
_spatial.cKDTree(coords[:s])
2021-04-29 12:26:40 +05:30
distances = tree.query(candidates)[0]
if distances.max() > distance: # require minimum separation
2021-11-02 22:31:32 +05:30
i = 0
coords[s] = candidates[distances.argmax()] # maximum separation to existing point cloud
2021-04-28 11:27:20 +05:30
s += 1
progress.update(s)
if i >= 100:
raise ValueError('seeding not possible')
2020-09-23 13:15:36 +05:30
2020-09-24 00:18:34 +05:30
return coords
2020-09-23 13:15:36 +05:30
def from_grid(grid,
selection: _IntCollection = None,
invert_selection: bool = False,
average: bool = False,
periodic: bool = True) -> _Tuple[_np.ndarray, _np.ndarray]:
2020-09-23 13:15:36 +05:30
"""
Create seeds from grid description.
2020-09-25 00:56:16 +05:30
2020-09-23 13:15:36 +05:30
Parameters
----------
2020-12-04 11:42:18 +05:30
grid : damask.Grid
Grid from which the material IDs are used as seeds.
selection : int or collection of int, optional
Material IDs to consider.
2022-01-13 03:43:38 +05:30
invert_selection : bool, optional
Consider all material IDs except those in selection. Defaults to False.
2022-01-13 03:43:38 +05:30
average : bool, optional
Seed corresponds to center of gravity of material ID cloud.
2022-01-13 03:43:38 +05:30
periodic : bool, optional
Center of gravity accounts for periodic boundaries.
2020-09-25 00:56:16 +05:30
Returns
-------
coords, materials : numpy.ndarray, shape (:,3); numpy.ndarray, shape (:)
2021-04-24 21:30:57 +05:30
Seed coordinates in 3D space, material IDs.
2020-09-23 13:15:36 +05:30
"""
2020-12-04 11:42:18 +05:30
material = grid.material.reshape((-1,1),order='F')
2022-03-09 20:05:36 +05:30
selection_ = _util.ensure_integer_list(selection)
mask = _np.full(grid.cells.prod(),True,dtype=bool) if selection_ is None else \
_np.isin(material,selection_,invert=invert_selection).flatten()
coords = _grid_filters.coordinates0_point(grid.cells,grid.size).reshape(-1,3,order='F')
2020-09-23 13:15:36 +05:30
if not average:
return (coords[mask],material[mask])
else:
materials = _np.unique(material[mask])
coords_ = _np.zeros((materials.size,3),dtype=float)
for i,mat in enumerate(materials):
2020-12-04 11:42:18 +05:30
pc = (2*_np.pi*coords[material[:,0]==mat,:]-grid.origin)/grid.size
coords_[i] = grid.origin + grid.size / 2 / _np.pi * (_np.pi +
_np.arctan2(-_np.average(_np.sin(pc),axis=0),
-_np.average(_np.cos(pc),axis=0))) \
if periodic else \
_np.average(coords[material[:,0]==mat,:],axis=0)
return (coords_,materials)