sisl.Grid
- class sisl.Grid(shape, bc=None, lattice=None, dtype=None, geometry=None)
Bases:
LatticeChild
,_Dispatchs
Real-space grid information with associated geometry.
This grid object handles cell vectors and divisions of said grid.
- Parameters:
shape (float or (3,) of int) – the shape of the grid. A
float
specifies the grid spacing in Angstrom, while a list of integers specifies the exact grid size.bc (list of int (3, 2) or (3, ), optional) – the boundary conditions for each of the cell’s planes. Default to periodic BC.
lattice (Lattice, optional) – the lattice that this grid represents.
lattice
has precedence if bothgeometry
andlattice
has been specified. Defaults to[1, 1, 1]
.dtype (numpy.dtype, optional) – the data-type of the grid, default to
numpy.float64
.geometry (Geometry, optional) – associated geometry with the grid. If
lattice
has not been passed the lattice will be taken from this geometry.
Examples
>>> grid1 = Grid(0.1, lattice=10) >>> grid2 = Grid(0.1, lattice=Lattice(10)) >>> grid3 = Grid(0.1, lattice=Lattice([10] * 3)) >>> grid1 == grid2 True >>> grid1 == grid3 True >>> grid = Grid(0.1, lattice=10, dtype=np.complex128) >>> grid == grid1 False
It is possible to provide a geometry and a different lattice to make a smaller (or bigger) lattice based on a geometry. This might be useful when creating wavefunctions or expanding densities to grids. Here we create a square grid based on a hexagonal graphene lattice. Expanding wavefunctions from this
geometry
will automatically convert to thelattice
size. >>> lattice = Lattice(10) # square lattice 10x10x10 Ang >>> geometry = geom.graphene() >>> grid = Grid(0.1, lattice=lattice, geometry=geometry)Methods
append
(other, axis)Appends other
Grid
to this grid along axisapply
(function_, *args, **kwargs)Applies a function to the grid and returns a new grid.
area
(ax0, ax1)Calculate the area spanned by the two axis ax0 and ax1 :rtype:
average
(axis[, weights])Average grid values along direction axis.
copy
([dtype])Copy the object, possibly changing the data-type
cross_section
(idx, axis)Takes a cross-section of the grid along axis axis
fill
(val)Fill the grid with this value
index
(coord[, axis])Find the grid index for a given coordinate (possibly only along a given lattice vector axis)
index2xyz
(index)Real-space coordinates of indices related to the grid
index_fold
(index[, unique])Converts indices from any placement to only exist in the "primary" grid
index_truncate
(index)Remove indices from outside the grid to only retain indices in the "primary" grid
interp
(shape[, order, mode])Interpolate grid values to a new grid of a different shape
isosurface
(level[, step_size])Calculates the isosurface for a given value.
mean
(axis[, weights])Average grid values along direction axis.
mgrid
(*slices)Return a list of indices corresponding to the slices
pyamg_boundary_condition
(A, b[, bc])Attach boundary conditions to the pyamg grid-matrix A with default boundary conditions as specified for this
Grid
pyamg_fix
(A, b, pyamg_indices, value)Fix values for the stencil to value.
pyamg_index
(index)Calculate pyamg matrix indices from a list of grid indices
pyamg_source
(b, pyamg_indices, value)Fix the source term to value.
read
(sile, *args, **kwargs)Reads grid from the
Sile
using read_gridremove
(idx, axis)Removes certain indices from a specified axis.
remove_part
(idx, axis, above)Removes parts of the grid via above/below designations.
sc_index
(*args, **kwargs)Call local
Lattice.sc_index
function :rtype:set_bc
(bc)set_boundary
(bc)set_geometry
(geometry[, also_lattice])Sets the
Geometry
for the grid.set_grid
(shape[, dtype])Create the internal grid of certain size.
set_lattice
(lattice)Overwrites the local lattice.
set_nsc
(*args, **kwargs)Set the number of super-cells in the
Lattice
objectset_sc
(lattice)Overwrites the local lattice.
set_supercell
(lattice)Overwrites the local lattice.
smooth
([r, method, mode])Make a smoother grid by applying a filter.
sub
(idx, axis)Retains certain indices from a specified axis.
sub_part
(idx, axis, above)Retains parts of the grid via above/below designations.
sum
(axis)Sum grid values along axis axis.
swapaxes
(a, b)Swap two axes in the grid (also swaps axes in the lattice)
tile
(reps, axis)Tile grid to create a bigger one
topyamg
([dtype])Create a pyamg stencil matrix to be used in pyamg
write
(sile, *args, **kwargs)Writes grid to the
Sile
using write_gridReturns the inherent
Lattice.cell
Voxel cell size
The data-type of the grid (in str)
Data-type used in grid
Volume of the grid voxel elements
Returns the inherent
Lattice.icell
Returns the inherent
Lattice.isc_off
Returns the inherent
Lattice.length
Returns the inherent
Lattice.n_s
Returns the inherent
Lattice.nsc
Returns the inherent
Lattice.origin
Returns the inherent
Lattice.rcell
[deprecated] Return the lattice object associated with the
Lattice
.Returns the inherent
Lattice.sc_off
Grid shape in \(x\), \(y\), \(z\) directions
Total number of elements in the grid
A dispatcher for classes, using __get__ it converts into ObjectDispatcher upon invocation from an object, or a TypeDispatcher when invoked from a class
Returns the inherent
Lattice.volume
- DIRICHLET = 3
- NEUMANN = 4
- OPEN = 5
- PERIODIC = 2
- apply(function_, *args, **kwargs)[source]
Applies a function to the grid and returns a new grid.
You can also apply a function that does not return a grid (maybe you want to do some measurement). In that case, you will get the result instead of a
Grid
.- Parameters:
function (str or function) – for a string the full module path to the function should be given. The function that will be called should have the grid as the first argument in its interface.
kwargs (args and) – arguments that go directly to the function call
- area(ax0, ax1)
Calculate the area spanned by the two axis ax0 and ax1 :rtype:
- average(axis, weights=None)[source]
Average grid values along direction axis.
- Parameters:
axis (int) – unit-cell direction to average across
weights (array_like, optional) – the weights for the individual axis elements, if boolean it corresponds to 0 and 1 for false/true.
See also
numpy.average
for details regarding the weights argument
- property cell: ndarray
Returns the inherent
Lattice.cell
- cross_section(idx, axis)[source]
Takes a cross-section of the grid along axis axis
Remark: This API entry might change to handle arbitrary cuts via rotation of the axis
- property dcell
Voxel cell size
- property dkind
The data-type of the grid (in str)
- property dtype
Data-type used in grid
- property dvolume
Volume of the grid voxel elements
- fill(val)[source]
Fill the grid with this value
- Parameters:
val (numpy.dtype) – all grid-points will have this value after execution
- property icell: ndarray
Returns the inherent
Lattice.icell
- index(coord, axis=None)[source]
Find the grid index for a given coordinate (possibly only along a given lattice vector axis)
- Parameters:
coord ((:, 3) or float or Shape) – the coordinate of the axis. If a float is passed axis is also required in which case it corresponds to the length along the lattice vector corresponding to axis. If a Shape a list of coordinates that fits the voxel positions are returned (all internal points also).
axis (int, optional) – the axis direction of the index, or for all axes if none.
- index2xyz(index)[source]
Real-space coordinates of indices related to the grid
- Parameters:
index (array_like) – indices for grid-positions
- Returns:
coordinates of the indices with respect to this grid spacing
- Return type:
- index_fold(index, unique=True)[source]
Converts indices from any placement to only exist in the “primary” grid
Examples
>>> grid = Grid([10, 10, 10]) >>> assert np.all(grid.index_fold([-1, -1, -1]) == 9)
- Parameters:
index (array_like) – indices for grid-positions
unique (bool, optional) – if true the returned indices are made unique after having folded the index points
- Returns:
all indices are then within the shape of the grid
- Return type:
See also
index_truncate
truncate indices by removing indices outside the primary cell
- index_truncate(index)[source]
Remove indices from outside the grid to only retain indices in the “primary” grid
Examples
>>> grid = Grid([10, 10, 10]) >>> assert len(grid.index_truncate([-1, -1, -1])) == 0
- Parameters:
index (array_like) – indices for grid-positions
- Returns:
all indices are then within the shape of the grid (others have been removed
- Return type:
See also
index_fold
fold indices into the primary cell
- interp(shape, order=1, mode='wrap', **kwargs)[source]
Interpolate grid values to a new grid of a different shape
It uses the
scipy.ndimage.zoom
, which creates a finer or more spaced grid using spline interpolation.- Parameters:
shape (int, array_like of len 3) – the new shape of the grid.
order (int 0-5, optional) – the order of the spline interpolation. 1 means linear, 2 quadratic, etc…
mode ({'wrap', 'mirror', 'constant', 'reflect', 'nearest'}) – determines how to compute the borders of the grid. The default is
'wrap'
, which accounts for periodic conditions.**kwargs – optional arguments passed to the interpolation algorithm The interpolation routine is
scipy.ndimage.zoom
See also
scipy.ndimage.zoom
method used for interpolation
- property isc_off: ndarray
Returns the inherent
Lattice.isc_off
- isosurface(level, step_size=1, **kwargs)[source]
Calculates the isosurface for a given value.
It uses
skimage.measure.marching_cubes
, so you need to have scikit-image installed.- Parameters:
level (float) – contour value to search for isosurfaces in the grid. If not given or None, the average of the min and max of the grid is used.
step_size (int, optional) – step size in voxels. Larger steps yield faster but coarser results. The result will always be topologically correct though.
**kwargs – optional arguments passed directly to
skimage.measure.marching_cubes
for the calculation of isosurfaces.
- Returns:
numpy array of shape (V, 3) – Verts. Spatial coordinates for V unique mesh vertices.
numpy array of shape (n_faces, 3) – Faces. Define triangular faces via referencing vertex indices from verts. This algorithm specifically outputs triangles, so each face has exactly three indices.
numpy array of shape (V, 3) – Normals. The normal direction at each vertex, as calculated from the data.
numpy array of shape (V, 3) – Values. Gives a measure for the maximum value of the data in the local region near each vertex. This can be used by visualization tools to apply a colormap to the mesh.
See also
skimage.measure.marching_cubes
method used to calculate the isosurface.
- property length: float
Returns the inherent
Lattice.length
- mean(axis, weights=None)
Average grid values along direction axis.
- Parameters:
axis (int) – unit-cell direction to average across
weights (array_like, optional) – the weights for the individual axis elements, if boolean it corresponds to 0 and 1 for false/true.
See also
numpy.average
for details regarding the weights argument
- classmethod mgrid(*slices)[source]
Return a list of indices corresponding to the slices
The returned values are equivalent to
numpy.mgrid
but they are returned in a (:, 3) array.
- property n_s: int
Returns the inherent
Lattice.n_s
- new = <TypeDispatcher{obj=<class 'sisl.Grid'>}>
- property nsc: ndarray
Returns the inherent
Lattice.nsc
- property origin: ndarray
Returns the inherent
Lattice.origin
- plot
Handles all plotting possibilities for a class
- pyamg_boundary_condition(A, b, bc=None)[source]
Attach boundary conditions to the pyamg grid-matrix A with default boundary conditions as specified for this
Grid
- Parameters:
A (scipy.sparse.csr_matrix) – sparse matrix describing the grid
b (numpy.ndarray) – a vector containing RHS of \(A x = b\) for the solution of the grid stencil
bc (list of BC, optional) – the specified boundary conditions. Default to the grid’s boundary conditions, else bc must be a list of elements with elements corresponding to
Grid.PERIODIC
/Grid.NEUMANN
…
- pyamg_fix(A, b, pyamg_indices, value)[source]
Fix values for the stencil to value.
- Parameters:
A (
csr_matrix
/csc_matrix
) – sparse matrix describing the LHS for the linear system of equationsb (numpy.ndarray) – a vector containing RHS of \(A x = b\) for the solution of the grid stencil
pyamg_indices (list of int) – the linear pyamg matrix indices where the value of the grid is fixed. I.e. the indices should correspond to returned quantities from pyamg_indices.
value (float) – the value of the grid to fix the value at
- pyamg_index(index)[source]
Calculate pyamg matrix indices from a list of grid indices
- Parameters:
index ((:, 3) of int) – a list of indices of the grid along each grid axis
- Returns:
linear indices for the matrix
- Return type:
See also
index
query indices from coordinates (directly passable to this method)
mgrid
Grid equivalent to
numpy.mgrid
. Grid.mgrid returns indices in shapes (:, 3), contrary to numpy’snumpy.mgrid
- Raises:
ValueError – if any of the passed indices are below 0 or above the number of elements per axis
- classmethod pyamg_source(b, pyamg_indices, value)[source]
Fix the source term to value.
- Parameters:
b (numpy.ndarray) – a vector containing RHS of \(A x = b\) for the solution of the grid stencil
pyamg_indices (list of int) – the linear pyamg matrix indices where the value of the grid is fixed. I.e. the indices should correspond to returned quantities from pyamg_indices.
- property rcell: ndarray
Returns the inherent
Lattice.rcell
- static read(sile, *args, **kwargs)[source]
Reads grid from the
Sile
using read_grid- Parameters:
sile (Sile, str or pathlib.Path) – a
Sile
object which will be used to read the grid if it is a string it will create a new sile usingget_sile
.* (args passed directly to
read_grid(,**)
) –
- remove(idx, axis)[source]
Removes certain indices from a specified axis.
Works exactly opposite to
sub
.- Parameters:
idx (array_like) – the indices of the grid axis axis to be removed
axis (int) – the axis segment from which we remove all indices idx
- remove_part(idx, axis, above)[source]
Removes parts of the grid via above/below designations.
Works exactly opposite to
sub_part
- sc_index(*args, **kwargs)
Call local
Lattice.sc_index
function :rtype:
- property sc_off: ndarray
Returns the inherent
Lattice.sc_off
- set_geometry(geometry, also_lattice=True)[source]
Sets the
Geometry
for the grid.Setting the
Geometry
for the grid is a possibility to attach atoms to the grid.It is not a necessary entity, so passing None is a viable option.
- set_lattice(lattice)
Overwrites the local lattice.
- set_nsc(*args, **kwargs)
Set the number of super-cells in the
Lattice
objectSee
set_nsc
for allowed parameters.See also
Lattice.set_nsc
the underlying called method
- set_sc(lattice)
Overwrites the local lattice.
- set_supercell(lattice)
Overwrites the local lattice.
- property shape
Grid shape in \(x\), \(y\), \(z\) directions
- property size
Total number of elements in the grid
- smooth(r=0.7, method='gaussian', mode='wrap', **kwargs)[source]
Make a smoother grid by applying a filter.
- Parameters:
r (float or array-like of len 3, optional) –
the radius of the filter in Angstrom for each axis. If the method is
"gaussian"
, this is the standard deviation!If a single float is provided, then the same distance will be used for all axes.
method ({'gaussian', 'uniform'}) – the type of filter to apply to smoothen the grid.
mode ({'wrap', 'mirror', 'constant', 'reflect', 'nearest'}) – determines how to compute the borders of the grid. The default is wrap, which accounts for periodic conditions.
See also
- sub(idx, axis)[source]
Retains certain indices from a specified axis.
Works exactly opposite to
remove
.- Parameters:
idx (array_like) – the indices of the grid axis axis to be retained
axis (int) – the axis segment from which we retain the indices idx
- sub_part(idx, axis, above)[source]
Retains parts of the grid via above/below designations.
Works exactly opposite to
remove_part
- sum(axis)[source]
Sum grid values along axis axis.
- Parameters:
axis (int) – unit-cell direction to sum across
- swapaxes(a, b)[source]
Swap two axes in the grid (also swaps axes in the lattice)
If
swapaxes(0,1)
it returns the 0 in the 1 values.
- tile(reps, axis)[source]
Tile grid to create a bigger one
The atomic indices for the base Geometry will be retained.
- Parameters:
:raises SislError : when the lattice is not commensurate with the geometry:
See also
Geometry.tile
equivalent method for Geometry class
- to
A dispatcher for classes, using __get__ it converts into ObjectDispatcher upon invocation from an object, or a TypeDispatcher when invoked from a class
This is a class-placeholder allowing a dispatcher to be a class attribute and converted into an ObjectDispatcher when invoked from an object.
If it is called on the class, it will return a TypeDispatcher.
This class should be an attribute of a class. It heavily relies on the __get__ special method.
- Parameters:
name (str) – name of the attribute in the class
dispatchs (dict, optional) – dictionary of dispatch methods
obj_getattr (callable, optional) – method with 2 arguments, an
obj
and theattr
which may be used to control how the attribute is called.instance_dispatcher (AbstractDispatcher, optional) – control how instance dispatchers are handled through __get__ method. This controls the dispatcher used if called from an instance.
type_dispatcher (AbstractDispatcher, optional) – control how class dispatchers are handled through __get__ method. This controls the dispatcher used if called from a class.
Examples
>>> class A: ... new = ClassDispatcher("new", obj_getattr=lambda obj, attr: getattr(obj.sub, attr))
The above defers any attributes to the contained A.sub attribute.
- topyamg(dtype=None)[source]
Create a pyamg stencil matrix to be used in pyamg
This allows retrieving the grid matrix equivalent of the real-space grid. Subsequently the returned matrix may be used in pyamg for solutions etc.
The pyamg suite is it-self a rather complicated code with many options. For details we refer to pyamg.
- Parameters:
dtype (numpy.dtype, optional) – data-type used for the sparse matrix, default to use the grid data-type
- Returns:
scipy.sparse.csr_matrix – the stencil for the pyamg solver
numpy.ndarray – RHS of the linear system of equations
Examples
This example proves the best method for a variety of cases in regards of the 3D Poisson problem:
>>> grid = Grid(0.01) >>> A, b = grid.topyamg() # automatically setups the current boundary conditions >>> # add terms etc. to A and/or b >>> import pyamg >>> from scipy.sparse.linalg import cg >>> ml = pyamg.aggregation.smoothed_aggregation_solver(A, max_levels=1000) >>> M = ml.aspreconditioner(cycle='W') # pre-conditioner >>> x, info = cg(A, b, tol=1e-12, M=M)
See also
pyamg_index
convert grid indices into the sparse matrix indices for
A
pyamg_fix
fixes stencil for indices and fixes the source for the RHS matrix (uses
pyamg_source
)pyamg_source
fix the RHS matrix
b
to a constant valuepyamg_boundary_condition
setup the sparse matrix
A
to given boundary conditions (called in this routine)
- property volume: float
Returns the inherent
Lattice.volume
- write(sile, *args, **kwargs)[source]
Writes grid to the
Sile
using write_grid- Parameters:
sile (Sile, str or pathlib.Path) – a
Sile
object which will be used to write the grid if it is a string it will create a new sile usingget_sile
- Return type: