Merge remote-tracking branch 'origin/development' into mpi-single-process-read
This commit is contained in:
commit
6dedaedea6
|
@ -564,8 +564,8 @@ class Geom:
|
|||
If given as floats, coordinates are addressed.
|
||||
center : int or float numpy.ndarray of shape (3)
|
||||
Center of the primitive. If given as integers, grid point
|
||||
locations (cell centers) are addressed.
|
||||
If given as floats, coordinates are addressed.
|
||||
coordinates (cell centers) are addressed.
|
||||
If given as floats, coordinates in space are addressed.
|
||||
exponent : numpy.ndarray of shape (3) or float
|
||||
Exponents for the three axes.
|
||||
0 gives octahedron (ǀxǀ^(2^0) + ǀyǀ^(2^0) + ǀzǀ^(2^0) < 1)
|
||||
|
@ -581,25 +581,27 @@ class Geom:
|
|||
Repeat primitive over boundaries. Defaults to True.
|
||||
|
||||
"""
|
||||
# normalized 'radius' and center
|
||||
r = np.array(dimension)/self.grid/2.0 if np.array(dimension).dtype in np.sctypes['int'] else \
|
||||
np.array(dimension)/self.size/2.0
|
||||
c = (np.array(center) + .5)/self.grid if np.array(center).dtype in np.sctypes['int'] else \
|
||||
(np.array(center) - self.origin)/self.size
|
||||
# radius and center
|
||||
r = np.array(dimension)/2.0*self.size/self.grid if np.array(dimension).dtype in np.sctypes['int'] else \
|
||||
np.array(dimension)/2.0
|
||||
c = (np.array(center) + .5)*self.size/self.grid if np.array(center).dtype in np.sctypes['int'] else \
|
||||
(np.array(center) - self.origin)
|
||||
|
||||
coords = grid_filters.cell_coord0(self.grid,np.ones(3)) \
|
||||
- ((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 = grid_filters.cell_coord0(self.grid,self.size,
|
||||
-(0.5*(self.size + (self.size/self.grid
|
||||
if np.array(center).dtype in np.sctypes['int'] else
|
||||
0)) if periodic else c))
|
||||
coords_rot = R.broadcast_to(tuple(self.grid))@coords
|
||||
|
||||
with np.errstate(all='ignore'):
|
||||
mask = np.sum(np.power(coords_rot/r,2.0**np.array(exponent)),axis=-1) > 1.0
|
||||
|
||||
if periodic: # translate back to center
|
||||
mask = np.roll(mask,((c-np.ones(3)*.5)*self.grid).astype(int),(0,1,2))
|
||||
mask = np.roll(mask,((c/self.size-0.5)*self.grid).round().astype(int),(0,1,2))
|
||||
|
||||
fill_ = np.full_like(self.material,np.nanmax(self.material)+1 if fill is None else fill)
|
||||
|
||||
return Geom(material = np.where(np.logical_not(mask) if inverse else mask, self.material,fill_),
|
||||
return Geom(material = np.where(np.logical_not(mask) if inverse else mask,
|
||||
self.material,
|
||||
np.nanmax(self.material)+1 if fill is None else fill),
|
||||
size = self.size,
|
||||
origin = self.origin,
|
||||
comments = self.comments+[util.execution_stamp('Geom','add_primitive')],
|
||||
|
|
|
@ -265,12 +265,12 @@ class TestGeom:
|
|||
@pytest.mark.parametrize('inverse',[True,False])
|
||||
@pytest.mark.parametrize('periodic',[True,False])
|
||||
def test_add_primitive_rotation(self,center,inverse,periodic):
|
||||
"""Rotation should not change result for sphere (except for discretization errors)."""
|
||||
g = np.array([32,32,32])
|
||||
"""Rotation should not change result for sphere."""
|
||||
g = np.random.randint(8,32,(3))
|
||||
s = np.random.random(3)+.5
|
||||
fill = np.random.randint(10)+2
|
||||
eu=np.array([np.random.randint(4),np.random.randint(2),np.random.randint(4)])*.5*np.pi
|
||||
G_1 = Geom(np.ones(g,'i'),[1.,1.,1.]).add_primitive(.3,center,1,fill,inverse=inverse,periodic=periodic)
|
||||
G_2 = Geom(np.ones(g,'i'),[1.,1.,1.]).add_primitive(.3,center,1,fill,Rotation.from_Eulers(eu),inverse,periodic=periodic)
|
||||
G_1 = Geom(np.ones(g,'i'),s).add_primitive(.3,center,1,fill,inverse=inverse,periodic=periodic)
|
||||
G_2 = Geom(np.ones(g,'i'),s).add_primitive(.3,center,1,fill,Rotation.from_random(),inverse,periodic=periodic)
|
||||
assert geom_equal(G_1,G_2)
|
||||
|
||||
|
||||
|
|
Loading…
Reference in New Issue