Mesh#
- class discretisedfield.Mesh(*, region=None, p1=None, p2=None, n=None, cell=None, bc='', subregions=None)#
Finite-difference mesh.
Mesh discretises the
discretisedfield.Region
, passed asregion
, using a regular finite-difference mesh. Since the region spans between two points \(\mathbf{p}_{1}\) and \(\mathbf{p}_{2}\), these points can be passed asp1
andp2
, instead of passingdiscretisedfield.Region
object. In this casediscretisedfield.Region
is created internally. Eitherregion
orp1
andp2
can be passed, not both. The region is discretised using a finite-difference cell, whose dimensions are defined withcell
. Alternatively, the domain can be discretised by passing the number of discretisation cellsn
in all three dimensions. Eithercell
orn
can be passed, not both.It is possible to define boundary conditions (bc) for the mesh by passing a string to
bc
.If it is necessary to define subregions in the mesh, a dictionary can be passed using
subregions
. More precisely, dictionary keys are strings (valid Python variable names), whereas values arediscretisedfield.Region
objects. It is necessary that subregions belong to the mesh region, are an aggregate of a discretisation cell, and are “aligned” with the mesh. If not,ValueError
is raised.In order to properly define a mesh, mesh region must be an aggregate of discretisation cells. Otherwise,
ValueError
is raised.- Parameters:
region (discretisedfield.Region, optional) – Cubic region to be discretised on a regular mesh. Either
region
orp1
andp2
should be defined, not both. Defaults toNone
.p2 (p1 /) – Diagonally-opposite region points, for example for three dimensions \(\mathbf{p} = (p_{x}, p_{y}, p_{z})\). Either
region
orp1
andp2
should be defined, not both. Defaults toNone
.cell (array_like, optional) – Discretisation cell size, for example for three dimensions \((d_{x}, d_{y}, d_{z})\). Either
cell
orn
should be defined, not both. Defaults toNone
.n (array_like, optional) – The number of discretisation cells, for example for three dimensions \((n_{x}, n_{y}, n_{z})\). Either
cell
orn
should be defined, not both. Defaults toNone
.bc (str, optional) – Periodic boundary conditions in geometrical directions. It is a string consisting of one or more characters representing the name of the direction(s) as present in
self.region.dims
, denoting the direction(s) along which the mesh is periodic. In the case of Neumann or Dirichlet boundary condition, string'neumann'
or'dirichlet'
is passed. Defaults to an empty string.subregions (dict, optional) – A dictionary defining subregions in the mesh. The keys of the dictionary are the region names (
str
) as valid Python variable names, whereas the values arediscretisedfield.Region
objects. Defaults to an empty dictionary.
- Raises:
ValueError – If mesh domain is not an aggregate of discretisation cells. Alternatively, if both
region
as well asp1
andp2
or bothcell
andn
are passed. Alternatively if one of the subregions is: (i) not in the mesh region, (ii) it is not an aggregate of discretisation cell, or (iii) it is not aligned with the mesh.
Examples
1. Defining a nano-sized thin film mesh by passing
region
andcell
parameters.>>> import discretisedfield as df ... >>> p1 = (-50e-9, -25e-9, 0) >>> p2 = (50e-9, 25e-9, 5e-9) >>> cell = (1e-9, 1e-9, 0.1e-9) >>> region = df.Region(p1=p1, p2=p2) >>> mesh = df.Mesh(region=region, cell=cell) >>> mesh Mesh(...)
2. Defining a nano-sized thin film mesh by passing
p1
,p2
andn
parameters.>>> n = (100, 50, 5) >>> mesh = df.Mesh(p1=p1, p2=p2, n=n) >>> mesh Mesh(...)
3. Defining a mesh with periodic boundary conditions in \(x\) and \(y\) directions.
>>> bc = 'xy' >>> region = df.Region(p1=p1, p2=p2) >>> mesh = df.Mesh(region=region, n=n, bc=bc) >>> mesh Mesh(...)
Defining a mesh with two subregions.
>>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> n = (10, 10, 10) >>> subregions = {'r1': df.Region(p1=(0, 0, 0), p2=(50, 100, 100)), ... 'r2': df.Region(p1=(50, 0, 0), p2=(100, 100, 100))} >>> mesh = df.Mesh(p1=p1, p2=p2, n=n, subregions=subregions) >>> mesh Mesh(...)
5. An attempt to define a mesh, whose region is not an aggregate of discretisation cells in the \(z\) direction.
>>> p1 = (-25, 3, 0) >>> p2 = (25, 6, 1) >>> cell = (5, 3, 0.4) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) Traceback (most recent call last): ... ValueError: ...
An attempt to define a mesh, whose subregion is not aligned.
>>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> cell = (10, 10, 10) >>> subregions = {'r1': df.Region(p1=(2, 0, 0), p2=(52, 100, 100))} >>> mesh = df.Mesh(p1=p1, p2=p2, subregions=subregions) Traceback (most recent call last): ... ValueError: ...
Methods
Extension of the
dir(self)
list.Relational operator
==
.Extracting the discretisation in a particular direction.
Extracts the mesh of a subregion.
Generator yielding coordinates of discretisation cells.
Number of discretisation cells in the mesh.
Return self|value.
Representation string.
Check if the mesh is close enough to the other based on a tolerance.
Axis selector.
Create a field whose values are the mesh coordinates.
Performs an N-dimensional discrete Fast Fourier Transform (FFT) on the mesh.
Performs an N-dimensional discrete inverse Fast Fourier Transform (iFFT) on the mesh.
Convert cell's index to its coordinate.
Check if meshes are aligned.
Line generator.
load_subregions
Load subregions from json file.
Mesh padding.
Convert point to the index of a cell which contains that point.
Slices of indices that correspond to cells contained in the region.
Rotate mesh by 90°.
save_subregions
Save subregions to json file.
Scale the underlying region and all subregions.
Select a part of the mesh.
Axis slider.
Translate the underlying region and all subregions.
Properties
Boundary condition for the mesh.
The cell size of the mesh.
Midpoints of the cells of the mesh along the spatial directions.
Discretisation cell volume.
Generator yielding indices of all mesh cells.
k3d
plot.matplotlib
plot.Number of cells along each dimension of the mesh.
pyvista
plot.Region on which the mesh is defined.
Subregions of the mesh.
Vertices of the cells of the mesh along the spatial directions.
- __dir__()#
Extension of the
dir(self)
list.For example in a three dimensional geometry with spatial dimensions
'x'
,'y'
, and'z'
, it adds'dx'
,'dy'
, and'dz'
.- Returns:
Avalilable attributes.
- Return type:
list
- __eq__(other)#
Relational operator
==
.Two meshes are considered to be equal if:
Regions of both meshes are equal.
Discretisation cell sizes are the same.
Boundary conditions
bc
andsubregions
are not considered to be necessary conditions for determining equality.- Parameters:
other (discretisedfield.Mesh) – Second operand.
- Returns:
True
if two meshes are equal andFalse
otherwise.- Return type:
bool
Examples
Check if meshes are equal.
>>> import discretisedfield as df ... >>> mesh1 = df.Mesh(p1=(0, 0, 0), p2=(5, 5, 5), cell=(1, 1, 1)) >>> mesh2 = df.Mesh(p1=(0, 0, 0), p2=(5, 5, 5), cell=(1, 1, 1)) >>> mesh3 = df.Mesh(p1=(1, 1, 1), p2=(5, 5, 5), cell=(2, 2, 2)) >>> mesh1 == mesh2 True >>> mesh1 != mesh2 False >>> mesh1 == mesh3 False >>> mesh1 != mesh3 True
- __getattr__(attr)#
Extracting the discretisation in a particular direction.
For example in a three dimensional geometry with spatial dimensions
'x'
,'y'
, and'z'
, if'dx'
,'dy'
, or'dz'
is accessed, the discretisation cell size in that direction is returned.- Parameters:
attr (str) – Discretisation direction (eg.
'dx'
,'dy'
, or'dz'
)- Returns:
Discretisation in a particular direction.
- Return type:
numbers.Real
Examples
Discretisation in the different directions.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> cell = (10, 25, 50) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) ... >>> mesh.dx 10.0 >>> mesh.dy 25.0 >>> mesh.dz 50.0
- __getitem__(item)#
Extracts the mesh of a subregion.
If subregions were defined by passing
subregions
dictionary when the mesh was created, this method returns a mesh defined on a subregion with keyitem
. Alternatively, adiscretisedfield.Region
object can be passed and a minimum-sized mesh containing it will be returned. The resulting mesh has the same discretisation cell as the original mesh. This method uses closed intervals, inclusive of endpoints i.e. [], for extracting the new mesh.- Parameters:
item (str, discretisedfield.Region) – The key of a subregion in
subregions
dictionary or a region object.- Returns:
Mesh of a subregion.
- Return type:
disretisedfield.Mesh
Example
Extract subregion mesh by passing a subregion key.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> cell = (10, 10, 10) >>> subregions = {'r1': df.Region(p1=(0, 0, 0), p2=(50, 100, 100)), ... 'r2': df.Region(p1=(50, 0, 0), p2=(100, 100, 100))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell, subregions=subregions) ... >>> len(mesh) # number of discretisation cells 1000 >>> mesh.region.pmin array([0, 0, 0]) >>> mesh.region.pmax array([100, 100, 100]) >>> submesh = mesh['r1'] >>> len(submesh) 500 >>> submesh.region.pmin array([0, 0, 0]) >>> submesh.region.pmax array([ 50, 100, 100])
Extracting a submesh on a “newly-defined” region.
>>> p1 = (-50e-9, -25e-9, 0) >>> p2 = (50e-9, 25e-9, 5e-9) >>> cell = (5e-9, 5e-9, 5e-9) >>> region = df.Region(p1=p1, p2=p2) >>> mesh = df.Mesh(region=region, cell=cell) ... >>> subregion = df.Region(p1=(0, 1e-9, 0), p2=(10e-9, 14e-9, 5e-9)) >>> submesh = mesh[subregion] >>> submesh.cell array([5.e-09, 5.e-09, 5.e-09]) >>> submesh.n array([2, 3, 1])
- __iter__()#
Generator yielding coordinates of discretisation cells.
The discretisation cell’s coordinate corresponds to its center point.
- Yields:
numpy.ndarray – For three dimensions, mesh cell’s center point \(\mathbf{p} = (p_{x}, p_{y}, p_{z})\).
Examples
Getting coordinates of all mesh cells.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (2, 2, 1) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(region=df.Region(p1=p1, p2=p2), cell=cell) >>> list(mesh) [array([0.5, 0.5, 0.5]), array([1.5, 0.5, 0.5]), array([0.5, 1.5, 0.5]),...]
See also
- __len__()#
Number of discretisation cells in the mesh.
It is computed by multiplying all elements of
n
:\[n_\text{total} = n_{x} n_{y} n_{z}.\]- Returns:
Total number of discretisation cells.
- Return type:
int
Examples
Getting the number of discretisation cells in a mesh.
>>> import discretisedfield as df ... >>> p1 = (0, 5, 0) >>> p2 = (5, 15, 2) >>> cell = (1, 0.1, 1) >>> mesh = df.Mesh(region=df.Region(p1=p1, p2=p2), cell=cell) >>> mesh.n array([ 5, 100, 2]) >>> len(mesh) 1000
- __or__(other)#
Return self|value.
- __repr__()#
Representation string.
Internally self._repr_html_() is called and all html tags are removed from this string.
- Returns:
Representation string.
- Return type:
str
Example
Getting representation string.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (2, 2, 1) >>> cell = (1, 1, 1) >>> bc = 'x' >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell, bc=bc) >>> mesh Mesh(Region(pmin=[0, 0, 0], pmax=[2, 2, 1], ...), n=[2, 2, 1], bc=x)
- allclose(other, rtol=None, atol=None)#
Check if the mesh is close enough to the other based on a tolerance.
This methods compares the two underlying regions using
Region.allclose
and the number of cellsn
of the two meshes. The value of relative tolerance (rtol
) and absolute tolerance (atol
) are passed on toRegion.allclose
for the comparison. If not provided default values ofRegion.allclose
are used.- Parameters:
other (discretisedfield.Mesh) – The other mesh used for comparison.
rtol (numbers.Real, optional) – Absolute tolerance. If
None
, the default value is the smallest edge length of the region multipled by theregion.tolerance_factor
.atol (numbers.Real, optional) – Relative tolerance. If
None
,region.tolerance_factor
is used.
- Returns:
True
if other mesh is close enough, otherwiseFalse
.- Return type:
bool
- Raises:
TypeError – If the
other
argument is not of typediscretisedfield.Mesh
or ifrtol
andatol
arguments are not of typenumbers.Real
.ValueError – If the dimensions of the mesh and the other mesh does not match.
Example
>>> p1 = (0, 0, 0) >>> p2 = (20e-9, 20e-9, 20e-9) >>> n = (10, 10, 10) >>> mesh1 = df.Mesh(p1=p1, p2=p2, n=n) ... >>> p1 = (0, 0, 0) >>> p2 = (20e-9 + 1.2e-12, 20e-9 + 1e-13, 20e-9 + 2e-12) >>> n = (10, 10, 10) >>> mesh2 = df.Mesh(p1=p1, p2=p2, n=n) ... >>> mesh1.allclose(mesh2, atol=1e-11) True >>> mesh1.allclose(mesh2, atol=1e-13) False
- axis_selector(*, widget='dropdown', description='axis')#
Axis selector.
For
widget='dropdown'
,ipywidgets.Dropdown
is returned, whereas forwidget='radiobuttons'
,ipywidgets.RadioButtons
is returned. Default widget description can be changed usingdescription
.- Parameters:
widget (str) – Type of widget to be returned. Defaults to
'dropdown'
.description (str) – Widget description to be showed. Defaults to
'axis'
.
- Returns:
Axis selection widget.
- Return type:
ipywidgets.Dropdown, ipywidgets.RadioButtons
Example
Get the
RadioButtons
slider.
>>> p1 = (0, 0, 0) >>> p2 = (10e-9, 10e-9, 10e-9) >>> n = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, n=n) ... >>> mesh.axis_selector(widget='radiobuttons') RadioButtons(...)
- coordinate_field()#
Create a field whose values are the mesh coordinates.
This method can be used to create a vector field with values equal to the coordinates of the cell midpoints. The result is equivalent to a field created with the following code:
mesh = df.Mesh(...) df.Field(mesh, dim=mesh.region.ndim, value=lambda point: point)
This method should be preferred over the manual creation with a callable because it provides much better performance.
- Returns:
Field with coordinates as values.
- Return type:
Examples
Create a coordinate field.
>>> import discretisedfield as df ... >>> mesh = df.Mesh(p1=(0, 0, 0), p2=(4, 2, 1), cell=(1, 1, 1)) >>> cfield = mesh.coordinate_field() >>> cfield Field(...)
Extract its value at position (0.5, 0.5, 0.5)
>>> cfield((0.5, 0.5, 0.5)) array([0.5, 0.5, 0.5])
Compare with manually created coordinate field
>>> manually = df.Field(mesh, nvdim=3, value=lambda point: point) >>> cfield.allclose(manually) True
- fftn(rfft=False)#
Performs an N-dimensional discrete Fast Fourier Transform (FFT) on the mesh.
This method computes the FFT in an N-dimensional space. The FFT is a way to transform a spatial-domain into a frequency domain. Note that any information about subregions in the mesh is lost during this transformation.
- Parameters:
rfft (bool, optional) – Determines if a real FFT is to be performed (if True) or a complex FFT (if False). Defaults to False, i.e., a complex FFT is performed by default.
- Returns:
A mesh representing the Fourier transform of the original mesh. The returned mesh has dimensions labeled with frequency (k) and cells have coordinates that correspond to the correct frequencies in the frequency domain.
- Return type:
Examples
1. Create a mesh and perform a FFT. >>> import discretisedfield as df >>> mesh = df.Mesh(p1=0, p2=10, cell=2) >>> fft_mesh = mesh.fftn() >>> fft_mesh.n array([5]) >>> fft_mesh.cell array([0.1]) >>> fft_mesh.region.pmin array([-0.25]) >>> fft_mesh.region.pmax array([0.25])
2. Perform a real FFT. >>> fft_mesh = mesh.fftn(rfft=True) >>> fft_mesh.n array([3]) >>> fft_mesh.cell array([0.1]) >>> fft_mesh.region.pmin array([-0.05]) >>> fft_mesh.region.pmax array([0.25])
3. Create a 2D mesh and perform a FFT. This demonstrates how the function works with higher dimensional meshes. >>> mesh = df.Mesh(p1=(0, 0), p2=(10, 10), cell=(2, 2)) >>> fft_mesh = mesh.fftn() >>> fft_mesh.n array([5, 5]) >>> fft_mesh.cell array([0.1, 0.1]) >>> fft_mesh.region.pmin array([-0.25, -0.25]) >>> fft_mesh.region.pmax array([0.25, 0.25])
- ifftn(rfft=False, shape=None)#
Performs an N-dimensional discrete inverse Fast Fourier Transform (iFFT) on the mesh.
This function calculates the iFFT in an N-dimensional space. The iFFT is a method to convert a frequency-domain signal into a spatial-domain signal. If ‘rfft’ is set to True and ‘shape’ is None, the original mesh shape is assumed to be even in the last dimension.
Please note that during Fourier transformations, the original position information is lost, causing the inverse Fourier transform to be centered at the origin. This can be rectified by mesh.translate to translate the mesh back to the desired position.
- Parameters:
rfft (bool, optional) – If set to True, a real FFT is performed. If False, a complex FFT is performed. Defaults to False.
shape ((tuple, np.ndarray, list), optional) – Specifies the shape of the original mesh. Defaults to None, which means the shape of the original mesh is used.
- Returns:
A mesh representing the inverse Fourier transform of the mesh.
- Return type:
Examples
1. Create a mesh and perform an iFFT. >>> import discretisedfield as df >>> mesh = df.Mesh(p1=0, p2=10, cell=2) >>> ifft_mesh = mesh.fftn().ifftn() >>> ifft_mesh.n array([5]) >>> ifft_mesh.cell array([2.]) >>> ifft_mesh.region.pmin array([-5.]) >>> ifft_mesh.region.pmax array([5.])
2. Perform a real iFFT. >>> ifft_mesh = mesh.fftn(rfft=True).ifftn(rfft=True, shape=mesh.n) >>> ifft_mesh.n array([5]) >>> ifft_mesh.cell array([2.]) >>> ifft_mesh.region.pmin array([-5.]) >>> ifft_mesh.region.pmax array([5.])
3. Perform a 2D iFFT. >>> mesh = df.Mesh(p1=(0, 0), p2=(10, 10), cell=(2, 2)) >>> ifft_mesh = mesh.fftn().ifftn() >>> ifft_mesh.n array([5, 5]) >>> ifft_mesh.cell array([2., 2.]) >>> ifft_mesh.region.pmin array([-5., -5.]) >>> ifft_mesh.region.pmax array([5., 5.])
4. Perform a real 2D iFFT. >>> ifft_mesh = mesh.fftn(rfft=True).ifftn(rfft=True, shape=mesh.n) >>> ifft_mesh.n array([5, 5]) >>> ifft_mesh.cell array([2., 2.]) >>> ifft_mesh.region.pmin array([-5., -5.]) >>> ifft_mesh.region.pmax array([5., 5.])
- index2point(index, /)#
Convert cell’s index to its coordinate.
- Parameters:
index (array_like) – For three dimensions, the cell’s index \((i_{x}, i_{y}, i_{z})\).
- Returns:
For three dimensions, the cell’s coordinate \(\mathbf{p} = (p_{x}, p_{y}, p_{z})\).
- Return type:
numpy.ndarray
- Raises:
ValueError – If
index
is out of range.
Examples
Converting cell’s index to its center point coordinate.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (2, 2, 1) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) >>> mesh.index2point((0, 0, 0)) array([0.5, 0.5, 0.5]) >>> mesh.index2point((0, 1, 0)) array([0.5, 1.5, 0.5])
See also
- is_aligned(other, tolerance=1e-12)#
Check if meshes are aligned.
Two meshes are considered to be aligned if and only if:
They have same discretisation cell size.
They have common cell coordinates.
for a given tolerance value.
- Parameters:
other (discretisedfield.Mesh) – Other mesh to be checked if it is aligned with self.
tolerance (int, float, optional) – The allowed extent of misalignment for discretisation cells and cell coordinates.
- Returns:
True
if meshes are aligned,False
otherwise.- Return type:
bool
- Raises:
TypeError – If
other
argument is not of typediscretisedfield.Mesh
or iftolerance
argument is not of typefloat
orint
.
Examples
Check if two meshes are aligned.
>>> import discretisedfield as df ... >>> p1 = (-50e-9, -25e-9, 0) >>> p2 = (50e-9, 25e-9, 5e-9) >>> cell = (5e-9, 5e-9, 5e-9) >>> region1 = df.Region(p1=p1, p2=p2) >>> mesh1 = df.Mesh(region=region1, cell=cell) ... >>> p1 = (-45e-9, -20e-9, 0) >>> p2 = (10e-9, 20e-9, 5e-9) >>> cell = (5e-9, 5e-9, 5e-9) >>> region2 = df.Region(p1=p1, p2=p2) >>> mesh2 = df.Mesh(region=region2, cell=cell) ... >>> p1 = (-42e-9, -20e-9, 0) >>> p2 = (13e-9, 20e-9, 5e-9) >>> cell = (5e-9, 5e-9, 5e-9) >>> region3 = df.Region(p1=p1, p2=p2) >>> mesh3 = df.Mesh(region=region3, cell=cell) ... >>> mesh1.is_aligned(mesh2) True >>> mesh1.is_aligned(mesh3) False >>> mesh1.is_aligned(mesh1) True >>> p_1 = (0, 0, 0) >>> p_2 = (0 + 1e-13, 0, 0) >>> p_3 = (0, 0, 0 + 1e-10) >>> p_4 = (20e-9, 20e-9, 20e-9) >>> p_5 = (20e-9 + 1e-13, 20e-9, 20e-9) >>> p_6 = (20e-9, 20e-9, 20e-9 + 1e-10) >>> cell = (5e-9, 5e-9, 5e-9) >>> mesh4 = df.Mesh(p1=p_1, p2=p_4, cell=cell) >>> mesh5 = df.Mesh(p1=p_2, p2=p_5, cell=cell) >>> mesh6 = df.Mesh(p1=p_3, p2=p_6, cell=cell) ... >>> mesh4.is_aligned(mesh5, 1e-12) True >>> mesh4.is_aligned(mesh6, 1e-11) False
- line(*, p1, p2, n)#
Line generator.
Given two points
p1
andp2
line is defined andn
points on that line are generated and yielded inn
iterations:\[\mathbf{r}_{i} = i\frac{\mathbf{p}_{2} - \mathbf{p}_{1}}{n-1}, \text{for}\, i = 0, ..., n-1\]- Parameters:
p2 (p1 /) – For three dimensions, points between which the line is defined \(\mathbf{p} = (p_{x}, p_{y}, p_{z})\).
n (int) – Number of points on the line.
- Yields:
tuple – \(\mathbf{r}_{i}\)
- Raises:
ValueError – If
p1
orp2
is outside the mesh region.
Examples
Creating line generator.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (2, 2, 2) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) ... >>> line = mesh.line(p1=(0, 0, 0), p2=(2, 0, 0), n=2) >>> list(line) [(0.0, 0.0, 0.0), (2.0, 0.0, 0.0)]
See also
plane()
- pad(pad_width)#
Mesh padding.
This method extends the mesh by adding (padding) discretisation cells in chosen direction(s). The way in which the mesh is going to be padded is defined by passing
pad_width
dictionary. The keys of the dictionary are the directions (axes), e.g.'x'
,'y'
, or'z'
, whereas the values are the tuples of length 2. The first integer in the tuple is the number of cells added in the negative direction, and the second integer is the number of cells added in the positive direction.- Parameters:
pad_width (dict) – The keys of the dictionary are the directions (axes), e.g.
'x'
,'y'
, or'z'
, whereas the values are the tuples of length 2. The first integer in the tuple is the number of cells added in the negative direction, and the second integer is the number of cells added in the positive direction.- Returns:
Padded (extended) mesh.
- Return type:
Examples
Padding a mesh in the x and y directions by 1 cell.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> cell = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) ... >>> mesh.region.edges array([100, 100, 100]) >>> padded_mesh = mesh.pad({'x': (1, 1), 'y': (1, 1), 'z': (0, 1)}) >>> padded_mesh.region.edges array([120., 120., 110.]) >>> padded_mesh.n array([12, 12, 11])
- point2index(point, /)#
Convert point to the index of a cell which contains that point.
This method uses half-open intervals for each cell, inclusive of the start point but exclusive of the endpoints. i.e. for each cell [). The exception to this is the very last cell contained in the region which has a closed interval i.e. [] and is inclusive of both the lower and upper bounds of the cell.
- Parameters:
point (array_like) – For three dimensions, point \(\mathbf{p} = (p_{x}, p_{y}, p_{z})\).
- Returns:
For three dimensions, the cell’s index \((i_{x}, i_{y}, i_{z})\).
- Return type:
tuple
- Raises:
ValueError – If
point
is outside the mesh.
Examples
Converting point to the cell’s index.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (2, 2, 1) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(region=df.Region(p1=p1, p2=p2), cell=cell) >>> mesh.point2index((0.2, 1.7, 0.3)) (0, 1, 0)
See also
- region2slices(region)#
Slices of indices that correspond to cells contained in the region.
- Parameters:
region (df.Region) – Region to convert to slices.
- Returns:
Tuple of slices of region indices.
- Return type:
tuple
Examples
1. Slices of a subregion >>> import discretisedfield as df … >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 1) >>> cell = (1, 1, 1) >>> subregions = {‘sr’: df.Region(p1=p1, p2=(10, 5, 1))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell, subregions=subregions) >>> mesh.region2slices(mesh.subregions[‘sr’]) (slice(0, 10, None), slice(0, 5, None), slice(0, 1, None))
- rotate90(ax1, ax2, k=1, reference_point=None, inplace=False)#
Rotate mesh by 90°.
Rotate the mesh
k
times by 90 degrees in the plane defined byax1
andax2
. The rotation direction is fromax1
toax2
, the two must be different.The rotate method does not rotate the string defining periodic boundary conditions, e.g. if a system has periodic boundary conditions in x and is rotated in the xy plane the new system will still have periodic boundary conditions in the new x direction, NOT in the new y direction. It is the user’s task to update the
bc
string after rotation if required.- Parameters:
ax1 (str) – Name of the first dimension.
ax2 (str) – Name of the second dimension.
k (int, optional) – Number of 90° rotations, defaults to 1.
reference_point (array_like, optional) – Point around which the mesh is rotated. If not provided the mesh.region’s centre point is used.
inplace (bool, optional) – If
True
, the rotation is applied in-place. Defaults toFalse
.
- Returns:
The rotated mesh object. Either a new object or a reference to the existing mesh for
inplace=True
.- Return type:
Examples
>>> import discretisedfield as df >>> import numpy as np >>> p1 = (0, 0, 0) >>> p2 = (10, 8, 6) >>> mesh = df.Mesh(p1=p1, p2=p2, n=(10, 4, 6)) >>> rotated = mesh.rotate90('x', 'y') >>> rotated.region.pmin array([ 1., -1., 0.]) >>> rotated.region.pmax array([9., 9., 6.]) >>> rotated.n array([ 4, 10, 6])
See also
- scale(factor, reference_point=None, inplace=False)#
Scale the underlying region and all subregions.
This method scales mesh.region and all subregions by a
factor
with respect to areference_point
. Iffactor
is a number the same scaling is applied along all dimensions. Iffactor
is array-like its length must matchregion.ndim
and different factors are applied along the different directions (based on their order). Ifreference_point
isNone
,mesh.region.center
is used as the reference point. A new object is created unlessinplace=True
is specified.Scaling the mesh also scales
mesh.cell
. The number of cellsmesh.n
stays constant.- Parameters:
factor (numbers.Real or array-like of numbers.Real) – Factor to scale the mesh.
reference_point (array_like, optional) – The position of the reference point is fixed when scaling the mesh. If not specified the mesh is scaled about its
mesh.region.center
.inplace (bool, optional) – If True, the mesh object is modified in-place. Defaults to False.
- Returns:
Resulting mesh.
- Return type:
- Raises:
ValueError, TypeError – If the operator cannot be applied.
Example
Scale a mesh without subregions.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1)) >>> res = mesh.scale(2) >>> res.region.pmin array([-5., -5., -5.]) >>> res.region.pmax array([15., 15., 15.])
Scale a mesh with subregions.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> sr = {'sub_reg': df.Region(p1=p1, p2=(5, 5, 5))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1), subregions=sr) >>> res = mesh.scale(2) >>> res.region.pmin array([-5., -5., -5.]) >>> res.region.pmax array([15., 15., 15.]) >>> res.subregions['sub_reg'].pmin array([-5., -5., -5.]) >>> res.subregions['sub_reg'].pmax array([5., 5., 5.])
Scale a mesh with subregions in place.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> sr = {'sub_reg': df.Region(p1=p1, p2=(5, 5, 5))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1), subregions=sr) >>> mesh.scale((2, 2, 5), inplace=True) Mesh(...) >>> mesh.region.pmin array([ -5., -5., -20.]) >>> mesh.region.pmax array([15., 15., 30.]) >>> mesh.subregions['sub_reg'].pmin array([ -5., -5., -20.]) >>> mesh.subregions['sub_reg'].pmax array([5., 5., 5.])
Scale with respect to the origin
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1)) >>> res = mesh.scale(2, reference_point=p1) >>> res.region.pmin array([0, 0, 0]) >>> res.region.pmax array([20, 20, 20])
See also
- sel(*args, **kwargs)#
Select a part of the mesh.
If one of the axis from
region.dims
is passed as a string, a mesh of a reduced dimension along the axis and perpendicular to it is extracted, intersecting the axis at its center. Alternatively, if a keyword (representing the axis) argument is passed with a real number value (e.g.x=1e-9
), a mesh of reduced dimensions intersects the axis at a point ‘nearest’ to the provided value is returned. If instead a tuple, list or a numpy array of length 2 is passed as a value containing two real numbers (e.g.x=(1e-9, 7e-9)
), a sub mesh is returned with minimum and maximum points along the selected axis, ‘nearest’ to the minimum and maximum of the selected values, respectively.- Parameters:
args – A string corresponding to the selection axis that belongs to
region.dims
.kwarg – A key corresponding to the selection axis that belongs to
region.dims
. The values are either anumbers.Real
or list, tuple, numpy array of length 2 containingnumbers.Real
which represents a point or a range of points to be selected from the mesh.
- Returns:
An extracted mesh.
- Return type:
Examples
Extracting the mesh at a specific point (
y=1
).
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (5, 5, 5) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) >>> mesh.region.ndim 3 >>> mesh.region.dims ('x', 'y', 'z') >>> plane_mesh = mesh.sel(y=1) >>> plane_mesh.region.ndim 2 >>> plane_mesh.region.dims ('x', 'z')
Extracting the xy-plane mesh at the mesh region center.
>>> plane_mesh = mesh.sel('z') >>> plane_mesh.region.ndim 2 >>> plane_mesh.region.dims ('x', 'y')
Specifying a range of points along axis
x
to be selected from mesh.
>>> selected_mesh = mesh.sel(x=(2, 4)) >>> selected_mesh.region.ndim 3 >>> selected_mesh.region.dims ('x', 'y', 'z')
- slider(axis, /, *, multiplier=None, description=None, **kwargs)#
Axis slider.
For
axis
, the name of a spatial dimension is passed. Based on that value,ipywidgets.SelectionSlider
is returned. Axis multiplier can be changed viamultiplier
.This method is based on
ipywidgets.SelectionSlider
, so any keyword argument accepted by it can be passed.- Parameters:
axis (str) – Axis for which the slider is returned (For eg.,
'x'
,'y'
, or'z'
).multiplier (numbers.Real, optional) – Axis multiplier. Defaults to
None
.
- Returns:
Axis slider.
- Return type:
ipywidgets.SelectionSlider
Example
Get the slider for the x-coordinate.
>>> p1 = (0, 0, 0) >>> p2 = (10e-9, 10e-9, 10e-9) >>> n = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, n=n) ... >>> mesh.slider('x') SelectionSlider(...)
- translate(vector, inplace=False)#
Translate the underlying region and all subregions.
This method translates mesh.region and all subregions by adding
vector
topmin
andpmax
. Thevector
must haveRegion.ndim
elements. A new object is created unlessinplace=True
is specified.- Parameters:
vector (array-like of numbers.Number) – Vector to translate the underlying region.
inplace (bool, optional) – If True, the Region objects are modified in-place. Defaults to False.
- Returns:
Resulting mesh.
- Return type:
- Raises:
ValueError, TypeError – If the operator cannot be applied.
Examples
Translate a mesh without subregions.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1)) >>> res = mesh.translate((2, -2, 5)) >>> res.region.pmin array([ 2, -2, 5]) >>> res.region.pmax array([12, 8, 15])
Translate a mesh with subregions.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> sr = {'sub_reg': df.Region(p1=p1, p2=(5, 5, 5))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1), subregions=sr) >>> res = mesh.translate((2, -2, 5)) >>> res.region.pmin array([ 2, -2, 5]) >>> res.region.pmax array([12, 8, 15]) >>> res.subregions['sub_reg'].pmin array([ 2, -2, 5]) >>> res.subregions['sub_reg'].pmax array([ 7, 3, 10])
Translate a mesh with subregions in place.
>>> import discretisedfield as df >>> p1 = (0, 0, 0) >>> p2 = (10, 10, 10) >>> sr = {'sub_reg': df.Region(p1=p1, p2=(5, 5, 5))} >>> mesh = df.Mesh(p1=p1, p2=p2, cell=(1, 1, 1), subregions=sr) >>> mesh.translate((2, -2, 5), inplace=True) Mesh(...) >>> mesh.region.pmin array([ 2, -2, 5]) >>> mesh.region.pmax array([12, 8, 15]) >>> mesh.subregions['sub_reg'].pmin array([ 2, -2, 5]) >>> mesh.subregions['sub_reg'].pmax array([ 7, 3, 10])
See also
- __hash__ = None#
- property bc#
Boundary condition for the mesh.
Periodic boundary conditions can be specified by passing a string containing one or more characters from
self.region.dims
(e.g.'x'
,'yz'
,'xyz'
for three dimensions). Neumann or Dirichlet boundary conditions are defined by passing'neumann'
or'dirichlet'
string. Neumann and Dirichlet boundary conditions are still experimental.- Returns:
A string representing periodic boundary condition along one or more axes, or Dirichlet or Neumann boundary condition. The string is empty if no boundary condition is defined.
- Return type:
str
- property cell#
The cell size of the mesh.
- Returns:
A numpy array representing discretisation size along respective axes.
- Return type:
numpy.ndarray
- property cells#
Midpoints of the cells of the mesh along the spatial directions.
This method returns a named tuple containing numpy arrays with midpoints of the cells along the spatial directions. Individual directions can be accessed from the tuple.
- Returns:
Namedtuple with elements corresponding to geometrical directions, the cell midpoints along the directions as numpy arrays.
- Return type:
collections.namedtuple
Examples
Getting midpoints along the
x
axis.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (10, 1, 1) >>> cell = (2, 1, 1) >>> mesh = df.Mesh(region=df.Region(p1=p1, p2=p2), cell=cell) ... >>> mesh.cells.x array([1., 3., 5., 7., 9.])
- property dV#
Discretisation cell volume.
- Returns:
Discretisation cell volume.
- Return type:
float
Examples
Discretisation cell volume.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> cell = (1, 2, 4) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) ... >>> mesh.dV 8.0
- property indices#
Generator yielding indices of all mesh cells.
- Yields:
tuple – For three dimensions, mesh cell indices \((i_{x}, i_{y}, i_{z})\).
Examples
Getting indices of all mesh cells.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (3, 2, 1) >>> cell = (1, 1, 1) >>> mesh = df.Mesh(p1=p1, p2=p2, cell=cell) >>> list(mesh.indices) [(0, 0, 0), (1, 0, 0), (2, 0, 0), (0, 1, 0), (1, 1, 0), (2, 1, 0)]
See also
- property k3d#
k3d
plot.If
plot
is not passed,k3d.Plot
object is created automatically. The color of the region and the discretisation cell can be specified usingcolor
length-2 tuple, where the first element is the colour of the region and the second element is the colour of the discretisation cell.It is often the case that the object size is either small (e.g. on a nanoscale) or very large (e.g. in units of kilometers). Accordingly,
multiplier
can be passed as \(10^{n}\), where \(n\) is a multiple of 3 (…, -6, -3, 0, 3, 6,…). According to that value, the axes will be scaled and appropriate units shown. For instance, ifmultiplier=1e-9
is passed, all axes will be divided by \(1\,\text{nm}\) and \(\text{nm}\) units will be used as axis labels. Ifmultiplier
is not passed, the best one is calculated internally.This method is based on
k3d.voxels
, so any keyword arguments accepted by it can be passed (e.g.wireframe
).- Parameters:
plot (k3d.Plot, optional) – Plot to which the plot is added. Defaults to
None
- plot is created internally.color ((2,) array_like) – Colour of the region and the discretisation cell. Defaults to the default color palette.
multiplier (numbers.Real, optional) – Axes multiplier. Defaults to
None
.
Examples
Visualising the mesh using
k3d
.
>>> p1 = (0, 0, 0) >>> p2 = (100, 100, 100) >>> n = (10, 10, 10) >>> mesh = df.Mesh(p1=p1, p2=p2, n=n) ... >>> mesh.k3d() Plot(...)
See also
- property mpl#
matplotlib
plot.If
ax
is not passed,matplotlib.axes.Axes
object is created automatically and the size of a figure can be specified usingfigsize
. The color of lines depicting the region and the discretisation cell can be specified usingcolor
length-2 tuple, where the first element is the colour of the region and the second element is the colour of the discretisation cell. The plot is saved in PDF-format iffilename
is passed.It is often the case that the object size is either small (e.g. on a nanoscale) or very large (e.g. in units of kilometers). Accordingly,
multiplier
can be passed as \(10^{n}\), where \(n\) is a multiple of 3 (…, -6, -3, 0, 3, 6,…). According to that value, the axes will be scaled and appropriate units shown. For instance, ifmultiplier=1e-9
is passed, all axes will be divided by \(1\,\text{nm}\) and \(\text{nm}\) units will be used as axis labels. Ifmultiplier
is not passed, the best one is calculated internally.This method is based on
matplotlib.pyplot.plot
, so any keyword arguments accepted by it can be passed (for instance,linewidth
,linestyle
, etc.).- Parameters:
ax (matplotlib.axes.Axes, optional) – Axes to which the plot is added. Defaults to
None
- axes are created internally.figsize ((2,) tuple, optional) – The size of a created figure if
ax
is not passed. Defaults toNone
.color ((2,) array_like) – A valid
matplotlib
color for lines depicting the region. Defaults to the default color palette.multiplier (numbers.Real, optional) – Axes multiplier. Defaults to
None
.box_aspect (str, array_like (3), optional) – Set the aspect-ratio of the plot. If set to ‘auto’ the aspect ratio is determined from the edge lengths of the region on which the mesh is defined. To set different aspect ratios a tuple can be passed. Defaults to
'auto'
.filename (str, optional) – If filename is passed, the plot is saved. Defaults to
None
.
Examples
Visualising the mesh using
matplotlib
.
>>> import discretisedfield as df ... >>> p1 = (-50e-9, -50e-9, 0) >>> p2 = (50e-9, 50e-9, 10e-9) >>> region = df.Region(p1=p1, p2=p2) >>> mesh = df.Mesh(region=region, n=(50, 50, 5)) ... >>> mesh.mpl()
See also
- property n#
Number of cells along each dimension of the mesh.
- Returns:
A numpy array representing number of discretisation cells along respective axes.
- Return type:
numpy.ndarray
- property pyvista#
pyvista
plot.
- property region#
Region on which the mesh is defined.
- Returns:
A region over which the regular mesh is defined.
- Return type:
- property subregions#
Subregions of the mesh.
When setting subregions all attributes of the individual regions (e.g. dims) apart from
pmin
andpmax
will be overwritten with the values frommesh.region
.- Returns:
A dictionary defining subregions in the mesh. The keys of the dictionary are the region names (
str
) as valid Python variable names, whereas the values arediscretisedfield.Region
objects.- Return type:
dict
- property vertices#
Vertices of the cells of the mesh along the spatial directions.
This method returns a named tuple containing numpy arrays with vertices of the cells along the spatial directions. Individual directions can be accessed from the tuple.
- Returns:
Namedtuple with elements corresponding to spatial directions, the cell vertices along the directions as numpy arrays.
- Return type:
collections.namedtuple
Examples
Getting vertices along the
x
axis.
>>> import discretisedfield as df ... >>> p1 = (0, 0, 0) >>> p2 = (10, 1, 1) >>> cell = (2, 1, 1) >>> mesh = df.Mesh(region=df.Region(p1=p1, p2=p2), cell=cell) ... >>> mesh.vertices.x array([ 0., 2., 4., 6., 8., 10.])