Skip to content

Commit

Permalink
plot geometry for dispersive materials without initializing structure…
Browse files Browse the repository at this point in the history
… object (NanoComp#1827)

* plot geometry without initializing structure class

* update docstrings

* rotate epsilon grid by 90 degrees to properly orient axes

* add support for dispersive ε

* update markdown pages from docstrings

* make frequency and resolution parameters of plot2D dictionary keys of eps_parameters

* reinstate frequency parameter of plot2D and add warning that it is deprecated

* fix order of frequency check
  • Loading branch information
oskooi authored and Mo Chen committed Feb 16, 2022
1 parent 9e3c8d5 commit d951cfd
Show file tree
Hide file tree
Showing 8 changed files with 270 additions and 185 deletions.
45 changes: 28 additions & 17 deletions doc/docs/Python_User_Interface.md
Original file line number Diff line number Diff line change
Expand Up @@ -676,16 +676,21 @@ frequency-independent part of $\mu$ (the $\omega\to\infty$ limit).
<div class="class_members" markdown="1">

```python
def get_epsilon_grid(self, xtics=None, ytics=None, ztics=None):
def get_epsilon_grid(self,
xtics=None,
ytics=None,
ztics=None,
frequency=0):
```

<div class="method_docstring" markdown="1">

Given three 1d NumPy arrays (`xtics`,`ytics`,`ztics`) which define the coordinates of a Cartesian
grid anywhere within the cell volume, compute the trace of the $\varepsilon$ tensor from the `geometry`
exactly at each grid point. (For [`MaterialGrid`](#materialgrid)s, the $\varepsilon$ at each grid
point is computed using bilinear interpolation from the nearest `MaterialGrid` points and possibly also
projected to form a level set.) Note that this is different from `get_epsilon_point` which computes
grid anywhere within the cell volume, compute the trace of the $\varepsilon(f)$ tensor at frequency
$f$ (in Meep units) from the `geometry` exactly at each grid point. `frequency` defaults to 0 which is
the instantaneous $\varepsilon$. (For [`MaterialGrid`](#materialgrid)s, the $\varepsilon$ at each
grid point is computed using bilinear interpolation from the nearest `MaterialGrid` points and possibly
also projected to form a level set.) Note that this is different from `get_epsilon_point` which computes
$\varepsilon$ by bilinearly interpolating from the nearest Yee grid points. This function is useful for
sampling the material geometry to any arbitrary resolution. The return value is a NumPy array with shape
equivalent to `numpy.meshgrid(xtics,ytics,ztics)`. Empty dimensions are collapsed.
Expand Down Expand Up @@ -2067,7 +2072,7 @@ including outside the cell and a `near2far` object, returns the computed
of fields $(E_x^1,E_y^1,E_z^1,H_x^1,H_y^1,H_z^1,E_x^2,E_y^2,E_z^2,H_x^2,H_y^2,H_z^2,...)$
in Cartesian coordinates and
$(E_r^1,E_\phi^1,E_z^1,H_r^1,H_\phi^1,H_z^1,E_r^2,E_\phi^2,E_z^2,H_r^2,H_\phi^2,H_z^2,...)$
in cylindrical coordinates for the frequencies 1,2,,`nfreq`.
in cylindrical coordinates for the frequencies 1,2,...,`nfreq`.

</div>

Expand Down Expand Up @@ -2609,7 +2614,7 @@ fr = mp.FluxRegion(volume=mp.GDSII_vol(fname, layer, zmin, zmax))

### Data Visualization

This module provides basic visualization functionality for the simulation domain. The spirit of the module is to provide functions that can be called with *no customization options whatsoever* and will do useful relevant things by default, but which can also be customized in cases where you *do* want to take the time to spruce up the output. The `Simulation` class provides the following methods:
This module provides basic visualization functionality for the simulation domain. The intent of the module is to provide functions that can be called with *no customization options whatsoever* and will do useful relevant things by default, but which can also be customized in cases where you *do* want to take the time to spruce up the output. The `Simulation` class provides the following methods:


<a id="Simulation.plot2D"></a>
Expand Down Expand Up @@ -2648,14 +2653,18 @@ sim.run(...)
field_func = lambda x: 20*np.log10(np.abs(x))
import matplotlib.pyplot as plt
sim.plot2D(fields=mp.Ez,
field_parameters={'alpha':0.8, 'cmap':'RdBu', 'interpolation':'none', 'post_process':field_func},
boundary_parameters={'hatch':'o', 'linewidth':1.5, 'facecolor':'y', 'edgecolor':'b', 'alpha':0.3})
field_parameters={'alpha':0.8, 'cmap':'RdBu', 'interpolation':'none', 'post_process':field_func},
boundary_parameters={'hatch':'o', 'linewidth':1.5, 'facecolor':'y', 'edgecolor':'b', 'alpha':0.3})
plt.show()
plt.savefig('sim_domain.png')
```
If you just want to quickly visualize the simulation domain without the fields (i.e., when
setting up your simulation), there is no need to invoke the `run` function prior to calling
`plot2D`. Just define the `Simulation` object followed by any DFT monitors and then
invoke `plot2D`.

Note: When running a [parallel simulation](Parallel_Meep.md), the `plot2D` function expects to be called
on all processes, but only generates a plot on the master process.
Note: When running a [parallel simulation](Parallel_Meep.md), the `plot2D` function expects
to be called on all processes, but only generates a plot on the master process.

**Parameters:**

Expand All @@ -2677,6 +2686,12 @@ on all processes, but only generates a plot on the master process.
- `alpha=1.0`: transparency of geometry
- `contour=False`: if `True`, plot a contour of the geometry rather than its image
- `contour_linewidth=1`: line width of the contour lines if `contour=True`
- `frequency=None`: for materials with a [frequency-dependent
permittivity](Materials.md#material-dispersion) $\varepsilon(f)$, specifies the
frequency $f$ (in Meep units) of the real part of the permittivity to use in the
plot. Defaults to the `frequency` parameter of the [Source](#source) object.
- `resolution=None`: the resolution of the $\varepsilon$ grid. Defaults to the
`resolution` of the `Simulation` object.
* `boundary_parameters`: a `dict` of optional plotting parameters that override
the default parameters for the boundary layers.
- `alpha=1.0`: transparency of boundary layers
Expand Down Expand Up @@ -2713,10 +2728,6 @@ on all processes, but only generates a plot on the master process.
- `alpha=0.6`: transparency of fields
- `post_process=np.real`: post processing function to apply to fields (must be
a function object)
* `frequency`: for materials with a [frequency-dependent
permittivity](Materials.md#material-dispersion) $\varepsilon(f)$, specifies the
frequency $f$ (in Meep units) of the real part of the permittivity to use in the
plot. Defaults to the `frequency` parameter of the [Source](#source) object.

</div>

Expand Down Expand Up @@ -4401,7 +4412,7 @@ def __init__(self,
medium2,
weights=None,
grid_type='U_DEFAULT',
do_averaging=False,
do_averaging=True,
beta=0,
eta=0.5,
damping=0):
Expand Down Expand Up @@ -6893,7 +6904,7 @@ A class used to record the fields during timestepping (i.e., a [`run`](#run-func
function). The object is initialized prior to timestepping by specifying the
simulation object and the field component. The object can then be passed to any
[step-function modifier](#step-function-modifiers). For example, one can record the
E<sub>z</sub> fields at every one time unit using:
$E_z$ fields at every one time unit using:

```py
animate = mp.Animate2D(sim,
Expand Down
2 changes: 1 addition & 1 deletion doc/docs/Python_User_Interface.md.in
Original file line number Diff line number Diff line change
Expand Up @@ -480,7 +480,7 @@ This feature is only available if Meep is built with [libGDSII](Build_From_Sourc

### Data Visualization

This module provides basic visualization functionality for the simulation domain. The spirit of the module is to provide functions that can be called with *no customization options whatsoever* and will do useful relevant things by default, but which can also be customized in cases where you *do* want to take the time to spruce up the output. The `Simulation` class provides the following methods:
This module provides basic visualization functionality for the simulation domain. The intent of the module is to provide functions that can be called with *no customization options whatsoever* and will do useful relevant things by default, but which can also be customized in cases where you *do* want to take the time to spruce up the output. The `Simulation` class provides the following methods:

@@ Simulation.plot2D @@
@@ Simulation.plot3D @@
Expand Down
12 changes: 7 additions & 5 deletions python/meep.i
Original file line number Diff line number Diff line change
Expand Up @@ -1095,12 +1095,12 @@ void _get_gradient(PyObject *grad, double scalegrad, PyObject *fields_a, PyObjec
$1 = (double *)array_data($input);
}

%typecheck(SWIG_TYPECHECK_POINTER, fragment="NumPy_Fragments") double* grid_vals {
%typecheck(SWIG_TYPECHECK_POINTER, fragment="NumPy_Fragments") std::complex<double>* grid_vals {
$1 = is_array($input);
}

%typemap(in, fragment="NumPy_Macros") double* grid_vals {
$1 = (double *)array_data($input);
%typemap(in, fragment="NumPy_Macros") std::complex<double>* grid_vals {
$1 = (std::complex<double> *)array_data($input);
}

// typemap for solve_cw:
Expand Down Expand Up @@ -2040,7 +2040,8 @@ void _get_epsilon_grid(geometric_object_list gobj_list,
int nx, double *xtics,
int ny, double *ytics,
int nz, double *ztics,
double *grid_vals) {
std::complex<double> *grid_vals,
double frequency) {
meep_geom::get_epsilon_grid(gobj_list,
mlist,
_default_material,
Expand All @@ -2051,7 +2052,8 @@ void _get_epsilon_grid(geometric_object_list gobj_list,
nx, xtics,
ny, ytics,
nz, ztics,
grid_vals);
grid_vals,
frequency);
}

%}
73 changes: 44 additions & 29 deletions python/simulation.py
Original file line number Diff line number Diff line change
Expand Up @@ -2225,18 +2225,19 @@ def get_mu_point(self, pt, frequency=0):
v3 = py_v3_to_vec(self.dimensions, pt, self.is_cylindrical)
return self.fields.get_mu(v3,frequency)

def get_epsilon_grid(self, xtics=None, ytics=None, ztics=None):
def get_epsilon_grid(self, xtics=None, ytics=None, ztics=None, frequency=0):
"""
Given three 1d NumPy arrays (`xtics`,`ytics`,`ztics`) which define the coordinates of a Cartesian
grid anywhere within the cell volume, compute the trace of the $\\varepsilon$ tensor from the `geometry`
exactly at each grid point. (For [`MaterialGrid`](#materialgrid)s, the $\\varepsilon$ at each grid
point is computed using bilinear interpolation from the nearest `MaterialGrid` points and possibly also
projected to form a level set.) Note that this is different from `get_epsilon_point` which computes
grid anywhere within the cell volume, compute the trace of the $\\varepsilon(f)$ tensor at frequency
$f$ (in Meep units) from the `geometry` exactly at each grid point. `frequency` defaults to 0 which is
the instantaneous $\\varepsilon$. (For [`MaterialGrid`](#materialgrid)s, the $\\varepsilon$ at each
grid point is computed using bilinear interpolation from the nearest `MaterialGrid` points and possibly
also projected to form a level set.) Note that this is different from `get_epsilon_point` which computes
$\\varepsilon$ by bilinearly interpolating from the nearest Yee grid points. This function is useful for
sampling the material geometry to any arbitrary resolution. The return value is a NumPy array with shape
equivalent to `numpy.meshgrid(xtics,ytics,ztics)`. Empty dimensions are collapsed.
"""
grid_vals = np.squeeze(np.empty((len(xtics), len(ytics), len(ztics))))
grid_vals = np.squeeze(np.empty((len(xtics), len(ytics), len(ztics)), dtype=np.complex128))
gv = self._create_grid_volume(False)
mp._get_epsilon_grid(self.geometry,
self.extra_materials,
Expand All @@ -2247,7 +2248,8 @@ def get_epsilon_grid(self, xtics=None, ytics=None, ztics=None):
len(xtics), xtics,
len(ytics), ytics,
len(ztics), ztics,
grid_vals)
grid_vals,
frequency)
return grid_vals

def get_filename_prefix(self):
Expand Down Expand Up @@ -2755,8 +2757,8 @@ def get_farfield(self, near2far, x):
(Fourier-transformed) "far" fields at `x` as list of length 6`nfreq`, consisting
of fields $(E_x^1,E_y^1,E_z^1,H_x^1,H_y^1,H_z^1,E_x^2,E_y^2,E_z^2,H_x^2,H_y^2,H_z^2,...)$
in Cartesian coordinates and
$(E_r^1,E_\phi^1,E_z^1,H_r^1,H_\phi^1,H_z^1,E_r^2,E_\phi^2,E_z^2,H_r^2,H_\phi^2,H_z^2,...)$
in cylindrical coordinates for the frequencies 1,2,,`nfreq`.
$(E_r^1,E_\\phi^1,E_z^1,H_r^1,H_\\phi^1,H_z^1,E_r^2,E_\\phi^2,E_z^2,H_r^2,H_\\phi^2,H_z^2,...)$
in cylindrical coordinates for the frequencies 1,2,...,`nfreq`.
"""
return mp._get_farfield(near2far.swigobj, py_v3_to_vec(self.dimensions, x, is_cylindrical=self.is_cylindrical))

Expand Down Expand Up @@ -4084,14 +4086,18 @@ def plot2D(self, ax=None, output_plane=None, fields=None, labels=False,
field_func = lambda x: 20*np.log10(np.abs(x))
import matplotlib.pyplot as plt
sim.plot2D(fields=mp.Ez,
field_parameters={'alpha':0.8, 'cmap':'RdBu', 'interpolation':'none', 'post_process':field_func},
boundary_parameters={'hatch':'o', 'linewidth':1.5, 'facecolor':'y', 'edgecolor':'b', 'alpha':0.3})
field_parameters={'alpha':0.8, 'cmap':'RdBu', 'interpolation':'none', 'post_process':field_func},
boundary_parameters={'hatch':'o', 'linewidth':1.5, 'facecolor':'y', 'edgecolor':'b', 'alpha':0.3})
plt.show()
plt.savefig('sim_domain.png')
```
If you just want to quickly visualize the simulation domain without the fields (i.e., when
setting up your simulation), there is no need to invoke the `run` function prior to calling
`plot2D`. Just define the `Simulation` object followed by any DFT monitors and then
invoke `plot2D`.
Note: When running a [parallel simulation](Parallel_Meep.md), the `plot2D` function expects to be called
on all processes, but only generates a plot on the master process.
Note: When running a [parallel simulation](Parallel_Meep.md), the `plot2D` function expects
to be called on all processes, but only generates a plot on the master process.
**Parameters:**
Expand All @@ -4113,6 +4119,12 @@ def plot2D(self, ax=None, output_plane=None, fields=None, labels=False,
- `alpha=1.0`: transparency of geometry
- `contour=False`: if `True`, plot a contour of the geometry rather than its image
- `contour_linewidth=1`: line width of the contour lines if `contour=True`
- `frequency=None`: for materials with a [frequency-dependent
permittivity](Materials.md#material-dispersion) $\\varepsilon(f)$, specifies the
frequency $f$ (in Meep units) of the real part of the permittivity to use in the
plot. Defaults to the `frequency` parameter of the [Source](#source) object.
- `resolution=None`: the resolution of the $\\varepsilon$ grid. Defaults to the
`resolution` of the `Simulation` object.
* `boundary_parameters`: a `dict` of optional plotting parameters that override
the default parameters for the boundary layers.
- `alpha=1.0`: transparency of boundary layers
Expand Down Expand Up @@ -4149,23 +4161,26 @@ def plot2D(self, ax=None, output_plane=None, fields=None, labels=False,
- `alpha=0.6`: transparency of fields
- `post_process=np.real`: post processing function to apply to fields (must be
a function object)
* `frequency`: for materials with a [frequency-dependent
permittivity](Materials.md#material-dispersion) $\\varepsilon(f)$, specifies the
frequency $f$ (in Meep units) of the real part of the permittivity to use in the
plot. Defaults to the `frequency` parameter of the [Source](#source) object.
"""
return vis.plot2D(self, ax=ax, output_plane=output_plane, fields=fields, labels=labels,
eps_parameters=eps_parameters, boundary_parameters=boundary_parameters,
source_parameters=source_parameters,monitor_parameters=monitor_parameters,
field_parameters=field_parameters, frequency=frequency,
plot_eps_flag=plot_eps_flag, plot_sources_flag=plot_sources_flag,
plot_monitors_flag=plot_monitors_flag, plot_boundaries_flag=plot_boundaries_flag,
**kwargs)


def plot_fields(self,**kwargs):
return vis.plot_fields(self,**kwargs)
return vis.plot2D(self,
ax=ax,
output_plane=output_plane,
fields=fields,
labels=labels,
eps_parameters=eps_parameters,
boundary_parameters=boundary_parameters,
source_parameters=source_parameters,
monitor_parameters=monitor_parameters,
field_parameters=field_parameters,
frequency=frequency,
plot_eps_flag=plot_eps_flag,
plot_sources_flag=plot_sources_flag,
plot_monitors_flag=plot_monitors_flag,
plot_boundaries_flag=plot_boundaries_flag,
**kwargs)

def plot_fields(self, **kwargs):
return vis.plot_fields(self, **kwargs)

def plot3D(self):
"""
Expand Down
2 changes: 1 addition & 1 deletion python/tests/test_visualization.py
Original file line number Diff line number Diff line change
Expand Up @@ -161,7 +161,7 @@ def test_animation_output(self):
sim = setup_sim() # generate 2D simulation

Animate = mp.Animate2D(sim,fields=mp.Ez, realtime=False, normalize=False) # Check without normalization
Animate_norm = mp.Animate2D(sim,mp.Ez,realtime=False,normalize=True) # Check with normalization
Animate_norm = mp.Animate2D(sim, mp.Ez, realtime=False, normalize=True) # Check with normalization

# test both animation objects during same run
sim.run(
Expand Down
Loading

0 comments on commit d951cfd

Please sign in to comment.