From 27d4512d8783629d70ffd5e7c7a6e01734ad24b7 Mon Sep 17 00:00:00 2001 From: Ardavan Oskooi Date: Thu, 31 Dec 2020 11:49:17 -0800 Subject: [PATCH] collapse empty dimensions of time-domain grid slices by default and remove snap argument from get_array_metadata (#1456) * fix for get_array_metadata related to collapsing of empty dimensions of grid slice * remove snap_empty_dimensions and collapse arrays in get_array_ * fix get_source_slice and get_array_metadata * update hard-coded validation binary files for cavity_arrayslice.py * remove collapse_empty_dimensions argument from tests/array_metadata.cpp * add array_to_all for parallel simulations and update docs * remove tests/array-slice-ll.cpp from make check suite because output_hdf5 only snaps to grid * add interpolation weights to slices and update hard-coded values of cavity_arrayslice.py test * re-instate snap argument for do_get_array_slice * fix bug in do_get_array_slice related to snapping * re-instate snap_empty_dimensions argument in get_array_slice_dimensions * include snap_empty_dimensions in meep.i * fixes * condense/simplify Co-authored-by: Steven G. Johnson --- doc/docs/Python_User_Interface.md | 50 ++-- doc/docs/Python_User_Interface.md.in | 2 +- python/simulation.py | 282 ++++++++++----------- python/tests/array_metadata.py | 7 +- python/tests/data/cavity_arrayslice_1d.npy | Bin 1088 -> 1136 bytes python/tests/data/cavity_arrayslice_2d.npy | Bin 38384 -> 38432 bytes python/tests/dispersive_eigenmode.py | 2 +- python/tests/faraday_rotation.py | 5 +- python/tests/simulation.py | 8 +- src/array_slice.cpp | 273 ++++++++++---------- src/meep.hpp | 32 +-- src/meepgeom.cpp | 3 +- tests/array-metadata.cpp | 17 +- tests/array-slice-ll-ref.h5 | Bin 42464 -> 42464 bytes tests/array-slice-ll.cpp | 8 +- 15 files changed, 336 insertions(+), 353 deletions(-) diff --git a/doc/docs/Python_User_Interface.md b/doc/docs/Python_User_Interface.md index 638f7146b..494f78e43 100644 --- a/doc/docs/Python_User_Interface.md +++ b/doc/docs/Python_User_Interface.md @@ -2457,8 +2457,8 @@ plt.show() plt.savefig('sim_domain.png') ``` -Note: When running a [parallel simulation](Parallel_Meep.md), the `plot2D` function must be called -from all parallel processes, but only produces a plot on the `meep.am_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:** @@ -3334,7 +3334,7 @@ See also [Field Functions](Field_Functions.md), and [Synchronizing the Magnetic #### Array Slices -The output functions described above write the data for the fields and materials for the entire cell to an HDF5 file. This is useful for post-processing as you can later read in the HDF5 file to obtain field/material data as a NumPy array. However, in some cases it is convenient to bypass the disk altogether to obtain the data *directly* in the form of a NumPy array without writing/reading HDF5 files. Additionally, you may want the field/material data on just a subregion (or slice) of the entire volume. This functionality is provided by the `get_array` method which takes as input a subregion of the cell and the field/material component. The method returns a NumPy array containing values of the field/material at the current simulation time. +The output functions described above write the data for the fields and materials for the entire cell to an HDF5 file. This is useful for post-processing large datasets which may not fit into memory as you can later read in the HDF5 file to obtain field/material data as a NumPy array. However, in some cases it is convenient to bypass the disk altogether to obtain the data *directly* in the form of a NumPy array without writing/reading HDF5 files. Additionally, you may want the field/material data on just a subregion (or slice) of the entire volume. This functionality is provided by the `get_array` method which takes as input a subregion of the cell and the field/material component. The method returns a NumPy array containing values of the field/material at the current simulation time. @@ -3350,7 +3350,7 @@ def get_array(self, cmplx=None, arr=None, frequency=0, - omega=0): + snap=False): ```
@@ -3370,7 +3370,7 @@ current simulation time. supplying a `Volume`. + `component`: field/material component (i.e., `mp.Ex`, `mp.Hy`, `mp.Sz`, - `mp.Dielectric`, etc). Defaults to `mp.Ez`. + `mp.Dielectric`, etc). Defaults to `None`. + `cmplx`: `boolean`; if `True`, return complex-valued data otherwise return real-valued data (default). @@ -3382,7 +3382,15 @@ current simulation time. fetching the same slice repeatedly at different times). + `frequency`: optional frequency point over which the average eigenvalue of the - dielectric and permeability tensors are evaluated (defaults to 0). + $\varepsilon$ and $\mu$ tensors are evaluated (defaults to 0). + ++ `snap`: By default, the elements of the grid slice are obtained using a bilinear + interpolation of the nearest Yee grid points. Empty dimensions of the grid slice + are "collapsed" into a single element. However, if `snap` is set to `True`, this + interpolation behavior is disabled and the grid slice is instead "snapped" + everywhere to the nearest grid point. (Empty slice dimensions are still of size + one.) This feature is mainly useful for comparing results with the + [`output_` routines](#output-functions) (e.g., `output_epsilon`, `output_efield_z`, etc.). For convenience, the following wrappers for `get_array` over the entire cell are available: `get_epsilon()`, `get_mu()`, `get_hpwr()`, `get_dpwr()`, @@ -3390,7 +3398,8 @@ available: `get_epsilon()`, `get_mu()`, `get_hpwr()`, `get_dpwr()`, `get_Xfield_z()`, `get_Xfield_r()`, `get_Xfield_p()` where `X` is one of `h`, `b`, `e`, `d`, or `s`. The routines `get_Xfield_*` all return an array type consistent with the fields (real or complex). The routines `get_epsilon()` and `get_mu()` -accept the optional `frequency` parameter (defaults to 0). +accept the optional argument `frequency` (defaults to 0) and all routines accept +`snap` (defaults to `False`). **Note on array-slice dimensions:** The routines `get_epsilon`, `get_Xfield_z`, etc. use as default `size=meep.Simulation.fields.total_volume()` which for @@ -3425,9 +3434,9 @@ Returns the Fourier-transformed fields as a NumPy array. + `dft_obj`: a `dft_flux`, `dft_force`, `dft_fields`, or `dft_near2far` object obtained from calling the appropriate `add` function (e.g., `mp.add_flux`). -+ `component`: a field component (e.g., `mp.Ez`) ++ `component`: a field component (e.g., `mp.Ez`). -+ `num_freq`: the index of the frequency: an integer in the range `0...nfreq-1`, ++ `num_freq`: the index of the frequency. An integer in the range `0...nfreq-1`, where `nfreq` is the number of frequencies stored in `dft_obj` as set by the `nfreq` parameter to `add_dft_fields`, `add_dft_flux`, etc. @@ -3449,8 +3458,6 @@ def get_array_metadata(self, center=None, size=None, dft_cell=None, - collapse=False, - snap=False, return_pw=False): ``` @@ -3462,7 +3469,7 @@ or `center`/`size`. In both cases, the return value is a tuple `(x,y,z,w)`, wher + `x,y,z` are 1d NumPy arrays storing the $x,y,z$ coordinates of the points in the grid slice -+ `w` is an array of the same dimensions as the array returned by ++ `w` is a NumPy array of the same dimensions as the array returned by `get_array`/`get_dft_array`, whose entries are the weights in a cubature rule for integrating over the spatial region (with the points in the cubature rule being just the grid points contained in the region). Thus, if $Q(\mathbf{x})$ is @@ -3472,15 +3479,23 @@ or `center`/`size`. In both cases, the return value is a tuple `(x,y,z,w)`, wher $$ \int_{\mathcal V} Q(\mathbf{x})d\mathbf{x} \approx \sum_{n} w_n Q_n.$$ This is a 1-, 2-, or 3-dimensional integral depending on the number of dimensions -in which $\mathcal{V}$ has zero extent. If the $\{Q_n\}$ samples are stored in an +in which $\mathcal{V}$ has zero extent. If the $Q_n$ samples are stored in an array `Q` of the same dimensions as `w`, then evaluating the sum on the RHS is just one line: `np.sum(w*Q).` A convenience parameter `dft_cell` is provided as an alternative to `vol` or -`center/size`; set `dft_cell` to a `dft_flux` or `dft_fields` object to define the -region covered by the array. If the `dft` argument is provided then all other -arguments (`vol`, `center`, and `size`) are ignored. If no arguments are provided, -then the entire cell is used. +`center`/`size`. Set `dft_cell` to a `dft_flux` or `dft_fields` object to define the +region covered by the array. If the `dft_cell` argument is provided then all other +arguments related to the spatial region (`vol`, `center`, and `size`) are ignored. +If no arguments are provided, then the entire cell is used. + +For empty dimensions of the grid slice `get_array_metadata` will collapse +the *two* elements corresponding to the nearest Yee grid points into a *single* +element using linear interpolation. + +If `return_pw=True`, the return value is a 2-tuple `(p,w)` where `p` (points) is a +list of `mp.Vector3`s with the same dimensions as `w` (weights). Otherwise, by +default the return value is a 4-tuple `(x,y,z,w)`.
@@ -3559,6 +3574,7 @@ def get_source(self, Return an array of complex values of the [source](#source) amplitude for `component` over the given `vol` or `center`/`size`. The array has the same dimensions as that returned by [`get_array`](#array-slices). +Not supported for [cylindrical coordinates](Python_Tutorials/Cylindrical_Coordinates.md). diff --git a/doc/docs/Python_User_Interface.md.in b/doc/docs/Python_User_Interface.md.in index 636e8bfe9..a90d8f16a 100644 --- a/doc/docs/Python_User_Interface.md.in +++ b/doc/docs/Python_User_Interface.md.in @@ -523,7 +523,7 @@ See also [Field Functions](Field_Functions.md), and [Synchronizing the Magnetic #### Array Slices -The output functions described above write the data for the fields and materials for the entire cell to an HDF5 file. This is useful for post-processing as you can later read in the HDF5 file to obtain field/material data as a NumPy array. However, in some cases it is convenient to bypass the disk altogether to obtain the data *directly* in the form of a NumPy array without writing/reading HDF5 files. Additionally, you may want the field/material data on just a subregion (or slice) of the entire volume. This functionality is provided by the `get_array` method which takes as input a subregion of the cell and the field/material component. The method returns a NumPy array containing values of the field/material at the current simulation time. +The output functions described above write the data for the fields and materials for the entire cell to an HDF5 file. This is useful for post-processing large datasets which may not fit into memory as you can later read in the HDF5 file to obtain field/material data as a NumPy array. However, in some cases it is convenient to bypass the disk altogether to obtain the data *directly* in the form of a NumPy array without writing/reading HDF5 files. Additionally, you may want the field/material data on just a subregion (or slice) of the entire volume. This functionality is provided by the `get_array` method which takes as input a subregion of the cell and the field/material component. The method returns a NumPy array containing values of the field/material at the current simulation time. @@ Simulation.get_array @@ @@ Simulation.get_dft_array @@ diff --git a/python/simulation.py b/python/simulation.py index 4f9781e7e..dc78fbd2e 100644 --- a/python/simulation.py +++ b/python/simulation.py @@ -3094,7 +3094,7 @@ def h5topng(self, rm_h5, option, *step_funcs): cmd = re.sub(r'\$EPS', self.last_eps_filename, opts) return convert_h5(rm_h5, cmd, *step_funcs) - def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None, arr=None, frequency=0, omega=0): + def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None, arr=None, frequency=0, snap=False): """ Takes as input a subregion of the cell and the field/material component. The method returns a NumPy array containing values of the field/material at the @@ -3111,7 +3111,7 @@ def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None supplying a `Volume`. + `component`: field/material component (i.e., `mp.Ex`, `mp.Hy`, `mp.Sz`, - `mp.Dielectric`, etc). Defaults to `mp.Ez`. + `mp.Dielectric`, etc). Defaults to `None`. + `cmplx`: `boolean`; if `True`, return complex-valued data otherwise return real-valued data (default). @@ -3123,7 +3123,15 @@ def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None fetching the same slice repeatedly at different times). + `frequency`: optional frequency point over which the average eigenvalue of the - dielectric and permeability tensors are evaluated (defaults to 0). + $\\varepsilon$ and $\\mu$ tensors are evaluated (defaults to 0). + + + `snap`: By default, the elements of the grid slice are obtained using a bilinear + interpolation of the nearest Yee grid points. Empty dimensions of the grid slice + are "collapsed" into a single element. However, if `snap` is set to `True`, this + interpolation behavior is disabled and the grid slice is instead "snapped" + everywhere to the nearest grid point. (Empty slice dimensions are still of size + one.) This feature is mainly useful for comparing results with the + [`output_` routines](#output-functions) (e.g., `output_epsilon`, `output_efield_z`, etc.). For convenience, the following wrappers for `get_array` over the entire cell are available: `get_epsilon()`, `get_mu()`, `get_hpwr()`, `get_dpwr()`, @@ -3131,7 +3139,8 @@ def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None `get_Xfield_z()`, `get_Xfield_r()`, `get_Xfield_p()` where `X` is one of `h`, `b`, `e`, `d`, or `s`. The routines `get_Xfield_*` all return an array type consistent with the fields (real or complex). The routines `get_epsilon()` and `get_mu()` - accept the optional `frequency` parameter (defaults to 0). + accept the optional argument `frequency` (defaults to 0) and all routines accept + `snap` (defaults to `False`). **Note on array-slice dimensions:** The routines `get_epsilon`, `get_Xfield_z`, etc. use as default `size=meep.Simulation.fields.total_volume()` which for @@ -3157,14 +3166,10 @@ def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None else: v = self._volume_from_kwargs(vol, center, size) - _, dirs = mp._get_array_slice_dimensions(self.fields, v, dim_sizes, False, True) + _, dirs = mp._get_array_slice_dimensions(self.fields, v, dim_sizes, not snap, snap) dims = [s for s in dim_sizes if s != 0] - if omega != 0: - frequency = omega - warnings.warn("get_array: omega has been deprecated; use frequency instead", RuntimeWarning) - if cmplx is None: cmplx = frequency != 0 or (component < mp.Dielectric and not self.fields.is_real) @@ -3183,9 +3188,9 @@ def get_array(self, component=None, vol=None, center=None, size=None, cmplx=None arr = np.zeros(dims, dtype=np.complex128 if cmplx else np.float64) if np.iscomplexobj(arr): - self.fields.get_complex_array_slice(v, component, arr, frequency) + self.fields.get_complex_array_slice(v, component, arr, frequency, snap) else: - self.fields.get_array_slice(v, component, arr, frequency) + self.fields.get_array_slice(v, component, arr, frequency, snap) return arr @@ -3198,9 +3203,9 @@ def get_dft_array(self, dft_obj, component, num_freq): + `dft_obj`: a `dft_flux`, `dft_force`, `dft_fields`, or `dft_near2far` object obtained from calling the appropriate `add` function (e.g., `mp.add_flux`). - + `component`: a field component (e.g., `mp.Ez`) + + `component`: a field component (e.g., `mp.Ez`). - + `num_freq`: the index of the frequency: an integer in the range `0...nfreq-1`, + + `num_freq`: the index of the frequency. An integer in the range `0...nfreq-1`, where `nfreq` is the number of frequencies stored in `dft_obj` as set by the `nfreq` parameter to `add_dft_fields`, `add_dft_flux`, etc. """ @@ -3225,23 +3230,21 @@ def get_source(self, component, vol=None, center=None, size=None): Return an array of complex values of the [source](#source) amplitude for `component` over the given `vol` or `center`/`size`. The array has the same dimensions as that returned by [`get_array`](#array-slices). + Not supported for [cylindrical coordinates](Python_Tutorials/Cylindrical_Coordinates.md). """ if vol is None and center is None and size is None: v = self.fields.total_volume() else: v = self._volume_from_kwargs(vol, center, size) dim_sizes = np.zeros(3, dtype=np.uintp) - mp._get_array_slice_dimensions(self.fields, v, dim_sizes, False, True) + mp._get_array_slice_dimensions(self.fields, v, dim_sizes, True, False) dims = [s for s in dim_sizes if s != 0] arr = np.zeros(dims, dtype=np.complex128) - self.fields.get_source_slice(v, component ,arr) + self.fields.get_source_slice(v, component, arr) return arr - # if return_pw, the return value is (points, weights) where points is a - # mp.Vector3-valued array of the same dimensions as weights. - # otherwise return value is 4-tuple (xtics, ytics, ztics, weights). def get_array_metadata(self, vol=None, center=None, size=None, dft_cell=None, - collapse=False, snap=False, return_pw=False): + return_pw=False): """ This routine provides geometric information useful for interpreting the arrays returned by `get_array` or `get_dft_array` for the spatial region defined by `vol` @@ -3249,7 +3252,7 @@ def get_array_metadata(self, vol=None, center=None, size=None, dft_cell=None, + `x,y,z` are 1d NumPy arrays storing the $x,y,z$ coordinates of the points in the grid slice - + `w` is an array of the same dimensions as the array returned by + + `w` is a NumPy array of the same dimensions as the array returned by `get_array`/`get_dft_array`, whose entries are the weights in a cubature rule for integrating over the spatial region (with the points in the cubature rule being just the grid points contained in the region). Thus, if $Q(\\mathbf{x})$ is @@ -3259,62 +3262,42 @@ def get_array_metadata(self, vol=None, center=None, size=None, dft_cell=None, $$ \\int_{\\mathcal V} Q(\\mathbf{x})d\\mathbf{x} \\approx \\sum_{n} w_n Q_n.$$ This is a 1-, 2-, or 3-dimensional integral depending on the number of dimensions - in which $\\mathcal{V}$ has zero extent. If the $\\{Q_n\\}$ samples are stored in an + in which $\\mathcal{V}$ has zero extent. If the $Q_n$ samples are stored in an array `Q` of the same dimensions as `w`, then evaluating the sum on the RHS is just one line: `np.sum(w*Q).` A convenience parameter `dft_cell` is provided as an alternative to `vol` or - `center/size`; set `dft_cell` to a `dft_flux` or `dft_fields` object to define the - region covered by the array. If the `dft` argument is provided then all other - arguments (`vol`, `center`, and `size`) are ignored. If no arguments are provided, - then the entire cell is used. + `center`/`size`. Set `dft_cell` to a `dft_flux` or `dft_fields` object to define the + region covered by the array. If the `dft_cell` argument is provided then all other + arguments related to the spatial region (`vol`, `center`, and `size`) are ignored. + If no arguments are provided, then the entire cell is used. + + For empty dimensions of the grid slice `get_array_metadata` will collapse + the *two* elements corresponding to the nearest Yee grid points into a *single* + element using linear interpolation. + + If `return_pw=True`, the return value is a 2-tuple `(p,w)` where `p` (points) is a + list of `mp.Vector3`s with the same dimensions as `w` (weights). Otherwise, by + default the return value is a 4-tuple `(x,y,z,w)`. """ if dft_cell: - vol, collapse = dft_cell.where, True + vol = dft_cell.where if vol is None and center is None and size is None: v = self.fields.total_volume() else: v = self._volume_from_kwargs(vol, center, size) - xyzw_vector = self.fields.get_array_metadata(v, collapse, snap) + xyzw_vector = self.fields.get_array_metadata(v) offset, tics = 0, [] for n in range(3): N = int(xyzw_vector[offset]) tics.append( xyzw_vector[offset+1:offset+1+N] ) offset += 1+N - wshape=[len(t) for t in tics if len(t)>1] - weights=np.reshape(xyzw_vector[offset:],wshape) + wshape = [len(t) for t in tics if len(t)>1] + weights = np.reshape(xyzw_vector[offset:], wshape) if return_pw: points=[ mp.Vector3(x,y,z) for x in tics[0] for y in tics[1] for z in tics[2] ] return points,weights - return tics + [weights] - - def get_dft_array_metadata(self, dft_cell=None, vol=None, center=None, size=None): - warnings.warn('get_dft_array_metadata is deprecated. Please use get_array_metadata instead', - DeprecationWarning) - return self.get_array_metadata(vol=dft_cell.where if dft_cell is not None else vol, - center=center, size=size, collapse=True) - - def get_array_slice_dimensions(self, component, vol=None, center=None, size=None): - """ - Computes the dimensions of a dft array for a particular `component` (`mp.Ez`, `mp.Ey`, etc.). - - Accepts either a volume object (`vol`), or a `center` and `size` `Vector3` pair. - - Returns a tuple containing the dimensions (`dim_sizes`), a `Vector3` object - corresponding to the minimum corner of the volume DFT volume object (`min_corner`), - and a `Vector3` object corresponding to the maximum corner (`max_corner`). - """ - if vol is None and center is None and size is None: - v = self.fields.total_volume() - else: - v = self._volume_from_kwargs(vol, center, size) - dim_sizes = np.zeros(3, dtype=np.uintp) - corners = [] - _,_ = mp._get_array_slice_dimensions(self.fields, v, dim_sizes, False, False, component,corners) - dim_sizes[dim_sizes==0] = 1 - min_corner = corners[0] - max_corner = corners[1] - return dim_sizes, min_corner, max_corner + return tuple(tics) + (weights,) def get_eigenmode_coefficients(self, flux, bands, eig_parity=mp.NO_PARITY, eig_vol=None, eig_resolution=0, eig_tolerance=1e-12, kpoint_func=None, direction=mp.AUTOMATIC): @@ -3675,153 +3658,150 @@ def output_times(self, fname): fname += '.csv' self.fields.output_times(fname) - def get_epsilon(self,frequency=0,omega=0): - if omega != 0: - frequency = omega - warnings.warn("get_epsilon: omega has been deprecated; use frequency instead", RuntimeWarning) - return self.get_array(component=mp.Dielectric,frequency=frequency) + def get_epsilon(self,frequency=0,snap=False): + return self.get_array(component=mp.Dielectric,frequency=frequency,snap=snap) - def get_mu(self): - return self.get_array(component=mp.Permeability) + def get_mu(self,frequency=0,snap=False): + return self.get_array(component=mp.Permeability,frequency=frequency,snap=snap) - def get_hpwr(self): - return self.get_array(component=mp.H_EnergyDensity) + def get_hpwr(self,snap=False): + return self.get_array(component=mp.H_EnergyDensity,snap=snap) - def get_dpwr(self): - return self.get_array(component=mp.D_EnergyDensity) + def get_dpwr(self,snap=False): + return self.get_array(component=mp.D_EnergyDensity,snap=snap) - def get_tot_pwr(self): - return self.get_array(component=mp.EnergyDensity) + def get_tot_pwr(self,snap=False): + return self.get_array(component=mp.EnergyDensity,snap=snap) - def get_hfield(self): + def get_hfield(self,snap=False): if self.is_cylindrical: - r = self.get_array(mp.Hr, cmplx=not self.fields.is_real) - p = self.get_array(mp.Hp, cmplx=not self.fields.is_real) + r = self.get_array(mp.Hr, cmplx=not self.fields.is_real, snap=snap) + p = self.get_array(mp.Hp, cmplx=not self.fields.is_real, snap=snap) return np.stack([r, p], axis=-1) else: - x = self.get_array(mp.Hx, cmplx=not self.fields.is_real) - y = self.get_array(mp.Hy, cmplx=not self.fields.is_real) - z = self.get_array(mp.Hz, cmplx=not self.fields.is_real) + x = self.get_array(mp.Hx, cmplx=not self.fields.is_real, snap=snap) + y = self.get_array(mp.Hy, cmplx=not self.fields.is_real, snap=snap) + z = self.get_array(mp.Hz, cmplx=not self.fields.is_real, snap=snap) return np.stack([x, y, z], axis=-1) - def get_hfield_x(self): - return self.get_array(mp.Hx, cmplx=not self.fields.is_real) + def get_hfield_x(self,snap=False): + return self.get_array(mp.Hx, cmplx=not self.fields.is_real, snap=snap) - def get_hfield_y(self): - return self.get_array(mp.Hy, cmplx=not self.fields.is_real) + def get_hfield_y(self,snap=False): + return self.get_array(mp.Hy, cmplx=not self.fields.is_real, snap=snap) - def get_hfield_z(self): - return self.get_array(mp.Hz, cmplx=not self.fields.is_real) + def get_hfield_z(self,snap=False): + return self.get_array(mp.Hz, cmplx=not self.fields.is_real, snap=snap) - def get_hfield_r(self): - return self.get_array(mp.Hr, cmplx=not self.fields.is_real) + def get_hfield_r(self,snap=False): + return self.get_array(mp.Hr, cmplx=not self.fields.is_real, snap=snap) - def get_hfield_p(self): - return self.get_array(mp.Hp, cmplx=not self.fields.is_real) + def get_hfield_p(self,snap=False): + return self.get_array(mp.Hp, cmplx=not self.fields.is_real, snap=snap) - def get_bfield(self): + def get_bfield(self,snap=False): if self.is_cylindrical: - r = self.get_array(mp.Br, cmplx=not self.fields.is_real) - p = self.get_array(mp.Bp, cmplx=not self.fields.is_real) + r = self.get_array(mp.Br, cmplx=not self.fields.is_real, snap=snap) + p = self.get_array(mp.Bp, cmplx=not self.fields.is_real, snap=snap) return np.stack([r, p], axis=-1) else: - x = self.get_array(mp.Bx, cmplx=not self.fields.is_real) - y = self.get_array(mp.By, cmplx=not self.fields.is_real) - z = self.get_array(mp.Bz, cmplx=not self.fields.is_real) + x = self.get_array(mp.Bx, cmplx=not self.fields.is_real, snap=snap) + y = self.get_array(mp.By, cmplx=not self.fields.is_real, snap=snap) + z = self.get_array(mp.Bz, cmplx=not self.fields.is_real, snap=snap) return np.stack([x, y, z], axis=-1) - def get_bfield_x(self): - return self.get_array(mp.Bx, cmplx=not self.fields.is_real) + def get_bfield_x(self,snap=False): + return self.get_array(mp.Bx, cmplx=not self.fields.is_real, snap=snap) - def get_bfield_y(self): - return self.get_array(mp.By, cmplx=not self.fields.is_real) + def get_bfield_y(self,snap=False): + return self.get_array(mp.By, cmplx=not self.fields.is_real, snap=snap) - def get_bfield_z(self): - return self.get_array(mp.Bz, cmplx=not self.fields.is_real) + def get_bfield_z(self,snap=False): + return self.get_array(mp.Bz, cmplx=not self.fields.is_real, snap=snap) - def get_bfield_r(self): - return self.get_array(mp.Br, cmplx=not self.fields.is_real) + def get_bfield_r(self,snap=False): + return self.get_array(mp.Br, cmplx=not self.fields.is_real, snap=snap) - def get_bfield_p(self): - return self.get_array(mp.Bp, cmplx=not self.fields.is_real) + def get_bfield_p(self,snap=False): + return self.get_array(mp.Bp, cmplx=not self.fields.is_real, snap=snap) - def get_efield(self): + def get_efield(self,snap=False): if self.is_cylindrical: - r = self.get_array(mp.Er, cmplx=not self.fields.is_real) - p = self.get_array(mp.Ep, cmplx=not self.fields.is_real) + r = self.get_array(mp.Er, cmplx=not self.fields.is_real, snap=snap) + p = self.get_array(mp.Ep, cmplx=not self.fields.is_real, snap=snap) return np.stack([r, p], axis=-1) else: - x = self.get_array(mp.Ex, cmplx=not self.fields.is_real) - y = self.get_array(mp.Ey, cmplx=not self.fields.is_real) - z = self.get_array(mp.Ez, cmplx=not self.fields.is_real) + x = self.get_array(mp.Ex, cmplx=not self.fields.is_real, snap=snap) + y = self.get_array(mp.Ey, cmplx=not self.fields.is_real, snap=snap) + z = self.get_array(mp.Ez, cmplx=not self.fields.is_real, snap=snap) return np.stack([x, y, z], axis=-1) - def get_efield_x(self): - return self.get_array(mp.Ex, cmplx=not self.fields.is_real) + def get_efield_x(self,snap=False): + return self.get_array(mp.Ex, cmplx=not self.fields.is_real, snap=snap) - def get_efield_y(self): - return self.get_array(mp.Ey, cmplx=not self.fields.is_real) + def get_efield_y(self,snap=False): + return self.get_array(mp.Ey, cmplx=not self.fields.is_real, snap=snap) - def get_efield_z(self): - return self.get_array(mp.Ez, cmplx=not self.fields.is_real) + def get_efield_z(self,snap=False): + return self.get_array(mp.Ez, cmplx=not self.fields.is_real, snap=snap) - def get_efield_r(self): - return self.get_array(mp.Er, cmplx=not self.fields.is_real) + def get_efield_r(self,snap=False): + return self.get_array(mp.Er, cmplx=not self.fields.is_real, snap=snap) - def get_efield_p(self): - return self.get_array(mp.Ep, cmplx=not self.fields.is_real) + def get_efield_p(self,snap=False): + return self.get_array(mp.Ep, cmplx=not self.fields.is_real, snap=snap) - def get_dfield(self): + def get_dfield(self,snap=False): if self.is_cylindrical: - r = self.get_array(mp.Dr, cmplx=not self.fields.is_real) - p = self.get_array(mp.Dp, cmplx=not self.fields.is_real) + r = self.get_array(mp.Dr, cmplx=not self.fields.is_real, snap=snap) + p = self.get_array(mp.Dp, cmplx=not self.fields.is_real, snap=snap) return np.stack([r, p], axis=-1) else: - x = self.get_array(mp.Dx, cmplx=not self.fields.is_real) - y = self.get_array(mp.Dy, cmplx=not self.fields.is_real) - z = self.get_array(mp.Dz, cmplx=not self.fields.is_real) + x = self.get_array(mp.Dx, cmplx=not self.fields.is_real, snap=snap) + y = self.get_array(mp.Dy, cmplx=not self.fields.is_real, snap=snap) + z = self.get_array(mp.Dz, cmplx=not self.fields.is_real, snap=snap) return np.stack([x, y, z], axis=-1) - def get_dfield_x(self): - return self.get_array(mp.Dx, cmplx=not self.fields.is_real) + def get_dfield_x(self,snap=False): + return self.get_array(mp.Dx, cmplx=not self.fields.is_real, snap=snap) - def get_dfield_y(self): - return self.get_array(mp.Dy, cmplx=not self.fields.is_real) + def get_dfield_y(self,snap=False): + return self.get_array(mp.Dy, cmplx=not self.fields.is_real, snap=snap) - def get_dfield_z(self): - return self.get_array(mp.Dz, cmplx=not self.fields.is_real) + def get_dfield_z(self,snap=False): + return self.get_array(mp.Dz, cmplx=not self.fields.is_real, snap=snap) - def get_dfield_r(self): - return self.get_array(mp.Dr, cmplx=not self.fields.is_real) + def get_dfield_r(self,snap=False): + return self.get_array(mp.Dr, cmplx=not self.fields.is_real, snap=snap) - def get_dfield_p(self): - return self.get_array(mp.Dp, cmplx=not self.fields.is_real) + def get_dfield_p(self,snap=False): + return self.get_array(mp.Dp, cmplx=not self.fields.is_real, snap=snap) - def get_sfield(self): + def get_sfield(self,snap=False): if self.is_cylindrical: - r = self.get_array(mp.Sr, cmplx=not self.fields.is_real) - p = self.get_array(mp.Sp, cmplx=not self.fields.is_real) + r = self.get_array(mp.Sr, cmplx=not self.fields.is_real, snap=snap) + p = self.get_array(mp.Sp, cmplx=not self.fields.is_real, snap=snap) return np.stack([r, p], axis=-1) else: - x = self.get_array(mp.Sx, cmplx=not self.fields.is_real) - y = self.get_array(mp.Sy, cmplx=not self.fields.is_real) - z = self.get_array(mp.Sz, cmplx=not self.fields.is_real) + x = self.get_array(mp.Sx, cmplx=not self.fields.is_real, snap=snap) + y = self.get_array(mp.Sy, cmplx=not self.fields.is_real, snap=snap) + z = self.get_array(mp.Sz, cmplx=not self.fields.is_real, snap=snap) return np.stack([x, y, z], axis=-1) - def get_sfield_x(self): - return self.get_array(mp.Sx, cmplx=not self.fields.is_real) + def get_sfield_x(self,snap=False): + return self.get_array(mp.Sx, cmplx=not self.fields.is_real, snap=snap) - def get_sfield_y(self): - return self.get_array(mp.Sy, cmplx=not self.fields.is_real) + def get_sfield_y(self,snap=False): + return self.get_array(mp.Sy, cmplx=not self.fields.is_real, snap=snap) - def get_sfield_z(self): - return self.get_array(mp.Sz, cmplx=not self.fields.is_real) + def get_sfield_z(self,snap=False): + return self.get_array(mp.Sz, cmplx=not self.fields.is_real, snap=snap) - def get_sfield_r(self): - return self.get_array(mp.Sr, cmplx=not self.fields.is_real) + def get_sfield_r(self,snap=False): + return self.get_array(mp.Sr, cmplx=not self.fields.is_real, snap=snap) - def get_sfield_p(self): - return self.get_array(mp.Sp, cmplx=not self.fields.is_real) + def get_sfield_p(self,snap=False): + return self.get_array(mp.Sp, cmplx=not self.fields.is_real, snap=snap) def plot2D(self, ax=None, output_plane=None, fields=None, labels=False, diff --git a/python/tests/array_metadata.py b/python/tests/array_metadata.py index 0fae2116a..83051292b 100644 --- a/python/tests/array_metadata.py +++ b/python/tests/array_metadata.py @@ -70,16 +70,13 @@ def vec_func(r): boundary_layers=pml_layers) dft_obj = sim.add_dft_fields([mp.Ez], fcen, 0, 1, where=nonpml_vol) - dft_obj_yee = sim.add_dft_fields([mp.Hy], fcen, 0, 1, where=nonpml_vol,yee_grid=True) sim.run(until_after_sources=100) - Hy_yee_dims,_,_ = sim.get_array_slice_dimensions(mp.Hy, dft_obj.where) - Ez = sim.get_dft_array(dft_obj, mp.Ez, 0) (X,Y,Z,W) = sim.get_array_metadata(dft_cell=dft_obj) - Eps = sim.get_array(vol=nonpml_vol,component=mp.Dielectric) + Eps = sim.get_array(vol=nonpml_vol, component=mp.Dielectric) EpsE2 = np.real(Eps*np.conj(Ez)*Ez) - xm, ym = np.meshgrid(X,Y) + xm, ym = np.meshgrid(X, Y) vec_func_sum = np.sum(W*(xm**2 + 2*ym**2)) pulse_modal_volume = np.sum(W*EpsE2)/np.max(EpsE2) * vec_func_sum diff --git a/python/tests/data/cavity_arrayslice_1d.npy b/python/tests/data/cavity_arrayslice_1d.npy index 056caeb073db63a09b165c1991f6d32a7ea0180c..9f23d858b612792b77dfdb5b2e8ae229afad59f4 100644 GIT binary patch literal 1136 zcmbW0|1;BZ9LJF@#I@++2d8y68+S`6#CE+u$L0ERb1WUHgqdpR%Nng8h|=ZTaT8)z zR71YpF}Yeoy?52JexYL=rLujPobzqtPR@1dy!H?D`Q>w;_kBK}_v`(9J|5Lk=Xuc= zY;B}AxvW%Swm`xPqqDea$60JTD=kZsBT3G@k|jwM{#%bsPR|x@uV-IN77MqZeF6f5 z*}iQ0ZTf#-PElrBgwp>x>cO1#sL}5zBC$WnxMoE0dUw8>Ln2Cb>}r67M6><<(&C3C zYJ97TDQ_lGFEe-TRU?V?;*|ELdnEc4mVbGJPa?^>@3&ZI5>@f-vy)Utlxhm^nBy4H zyQaC`&>jPtQBR)s*>6Az3hCFASEtdifqO2aVFFE$)-Mmv%8{r-?JX?U<8fIVrAs`6 z{p8ZpjD$seShx3B9M6n1GI@0#VHT{JR}_XvS#Vi}yt|NR!HngWH%CG&xJ1-H%h+SV zr-C~x5f(-iS@D_FfGd9f9;f?x7G&M4=W z@asO+#u#T2*VfD26?qF-r=Cn}NHE~{`D#WBQ;%kTo_H21TSbno0$rtW3puzX`OlZz z5zbsPuc(wdg}|8_xDx}1uE3=ya5@FtMgqs|@g1%mj60mi0rw{Gpa;B20#C}on{@DK e2)r5s&-%bStA{4=(iS{T0dK7y4}jNa!SlZwd?JSc literal 1088 zcmX|=e=L-79LN2*Gs05LkJco+b1W-oA@0`qnZ_oaaj{VO(cQ_e{5o~o$gf25<3bHL zqQm@XoTytU-`jMh^#`jXRO9YRml?+$o0(dgPd)!Ud%k;~=ll7-->+9)L<}b~&eSB- zB$vS#Wb;IfAU6g(0NYwgT z9#z>+qURpD8!sM{Xgc#sciSBjy$#C0yu~F^)~45&Xd4n$b1B(LaxD`62=39bv}m+V z_cCxmgS3j-Fi!^!ij#`xLar{L2@U&vTJsEY_|&L>x7dLARSLSGOpUdYE^1%qckI<5 zE>DYJ!Hn_!CogaeSeVYKcMdY*u_bA7aD)+8R5kP$bBy?yzGLV_fDso6M;2+;MjYbP5AMo9fV0^Zb_5;?XiPP4>>$WY-;&B|}Y zUa+v8?%Y9q*-=jeU5})^?bhW2Bd(DIPWskZDQzlqs2|Jil{+lCUF3qZa^Rjix1sj1 zk}2-2SBhPgfpG7G^)+~p0^ha6e@fuLL9h=C_RWJiZoynlFy}nXt%Gwg;avW3&VT2g zf*eXA7Z=Eh4!Jc!j`onNE#zzpx$B{aF6c!KJyk<*ub{`L(CaYt>;kAIHfyH6xC<-YpY8ufY zjC2wP!>A@@8DtnshN}Zsb@Pa zS>){2yH`Y}x}##m|M^ejng6?8Zx!IRBoX#8a z{-ifcHJa?55+q~YK6c|y++4wIceW53YfrMwV{cg}uOoeTN$)&MR>YIMmuG39e+o%{ zFn#tb(=2i_)H-QxS`K-y)hB)8lRT1gq-~=2=mIi*uW5C6ejzy#5S(4T_W+SD-5dB+ zRzeoO81!Y)q*CIHPfz@vg8LfZ3&nMw;=0>E+(*H0v<5d$!9aR)Ka4D%=wPJv*<$29 zn8{YmY!_x)ika`jN?gRsoW@FNjY+yN)!m0#-TcaCheH~RUuQJ9yfTd)%v%vyDC_%68`u*#~Ve&OQ4bNNT zZCd}RSoJ%kf?2L3@5VH$nV0IYj ztWYuV7)<0OCKiW@2Jo4Gn0ViHqs3B~U}+whi>0bf&iQBFN*`8ud|20$TPc*O>GYP2 z5`vwAnfoxi|2_>b_S|G~T=;KEw6hO~ZV16B%qt}$+rz9p*k{9tr~B!6D(T0eyM`*k z4#Ny+pY~0ns~fv9zW9BIwSr3fbI5OkN-!^&pDau-os_bS<)=%WPq}QQiUAyYAq3Nf zdDUm0S!~Zz_WQPu{*7-p(5g>3J0c#HphOGPI&K4@y;?Qm(*f3*ZPBXXEdlc@=cFo`C=~S0S3kGs%lMu`W zW|i%NryAt}Z1KvX!~xPQ+NjN;I3d{2Fgsra6ggW2GasusHfdxIRSx2iw-D?p%(w|B z&S*Vc&4SCn$<~;XM=f+Y6ea{Kf~iPrZr|CohFv?|rlq-jAC(Q}&}|`DGR$A2X6)6? z4`W8HeS_Sm6w(Ym4w(zV)HVon;k?Mua;x>M>1gqQ*SZHN8^WRcLa>!EU2~_X+hm8c z-@3fMxk-y@=THuH3BlrFZdsfq4UN=C?4j}z9C8tY z-QO+Dnb9lm9ntq)4|N+Gxrl%^NDg zW`+y1yX49dgO0hxA%3D}rr}bm9L1q0;VMCbiiLS*%n1|wE+1l)WwGd~PAat+&7moo zDnV>vK7S_j3J70A^xv8uJ=vH~WuJ3M@kS+x4Cc_tsGNphBFX)+J!6e$mCy_`4%G=k z&P52bBcih1B~U>woa;4mcuOT^<{XLcxvs4Lk&{&waGp8j;6~U}eNF}X%oI|CALp8xF!N$R~dNN9{pMMNem&rod#4#|aJJ79+ThmE)E-ogB|vwJ)qoTFq8 zhvI}_i7^)NgBQaLcs;r%lZ?<1AFap<&oQ^|3dhwvsF;x#W=!J#pD z&Ff(rV^7GzhSBKGp<-+p!(qP1*3*EkN9n<#TiANu!90h3YAH4`3r`Nou!&uT>5Cn& z7~7#t#-T&l4x?cHdt=S;;!rF$);-u*L$M3Xu|+d)4z0x&oeeVz`*{;K`%WJYL6BIO z({M7l;n9r2ad98Hi*Jc9fe=UIu50x@Vf;w7**oA08v7+^&ILEf*poAANA&V6r-Kt9Ew0OS_bo0 zx+rTsP}WK!Ih2XAR)VtD0X47$1+p}XL*6KmlVDb%nqGoZTd{#dktnsF!TeaX_uVUs z_Nr(OjT3_9qiEM+qWmvM`QN^gLuDxc4PahJui%14Li<||1)`C-0rLlR9rk~TmW0G| zsJ{@b1ZElfm(P!h=0>jI(0epD*)XS~6LLeFR1n7@TeL~TVAi0gno}nluEtFq(h!1W z!#sxWEE}zxGM+2pg)T+U!8ZxbA;Oy|d|Ak}7v?k(4gK*>3 z7>E%Dwm>ydMK#&7T7+6AAz~(%dqi4igqc3uEkdoDFp+8PzaiqjTCCJhhyx!X9;`%M zSdaMdBI3krh!<^qMR2+%i^q=&o}B6-kMvU|r}4z&ZxNizkYAP}IBi70nS{V|76Q;( z1fmNNkbZ^0^kaZJ8-FhZXNm89&#$xoLjv1`3-8B`a&Y5X1g8}km=ps&yh8-1_x=@} TredbsG4p7w#38Iq7FO!N!2ze8 delta 5770 zcmZ9Qd0fp|AIFJiETt)-Azcz0Tckq8$(AHaQW#XyqVbenJWaPD?JC_CQf|_wiI!8; zD5p(KLyLWzxnZoYOeQ^Kz0a9@KK_}%?(1_u@B90mbMEhV&iDJtEL3~-vs!=IvXx7n zR`>2@d#JyybU=sq#~!~6sjP6kmvwMN3M*Fo`n|0|G7Fz^R;fH7V|CSOLq1m?VzDEp zdB*s~v+fOZCsy2vW;V&!2o2dr5{*OdT1A~8eGfbS-Vs$w(*9RZeZ9|-s7DsF-x}7D z$^fhItp^*)XU#rwQ=T=En39`QJjb>Wiwwik?xWH+QsLt-E6un@99=VfU%1~O4qe~O z&zgFZIN@_WzNg@ptA5x<0&t-hxNrh)v;{Yg!9aRa#lW;N&>EQ-c+Wo*k@A^DVrGXi z(_GB_7*?VND^rb?YJ6P1JnVfEi@Utzk^V>-v%kLTNmq9QOZRm>YxpXbdDJ_sRhUIE zqtOpqe@F{ut#dBlIr3&JbFsK-7qBUfr2p534OcHGMoSeX(zp|~#93*wC^O{}DNnbR z9kFX6^A_t*|8(;zDP7_ebWL%C>ix3g-WdP|J) z{W)-{k7ZvS8<4dt5Y?^wsVyu@Qn%n3r=#E4IywW$h#F5=#6JP)BtRB?!Ta zVfrY{HifxFG1ED-l*8`)M5SMGXrB(M=NTW%wzS_3$Ab%;u^?~8@` zSUSx1abXnmI6hv}`om9@4CGL~5Nv9wFo};=OzPjUY|TN_`e%<0&_E3i-4%lEfq7%Z zIm7KVf!)q)u>8C*hRQWL6fFcBt}je8V{_dnJ7p|$nYKkgWddy(#37jwYz53Bn zoJ?V};w#(JI+CbTi$eiIuq!Y-yL_^oOw!nrqI-8k?T=6sZ4TKB!Ct@&om^3?`8b0O z%rF_M@hXG5>u^X{2$ls?dazBA+?vViK5Hz0=z5&S59ZKCAy^d5dqFcYbdF}R_syw) z&d)nRnJ$O=x2ObD+apW|yWIgNXJxbO=Jw^$mvd<65Drxe!PdcaNj2?nEz4o!ngW7) z^7E)<7>8yE!S=z7aMRJQipyo2qie3z8l0q#dK~f+f@Q&ccV4|`pl%-Xb*yL`sa-&& z-*U)J2=)qQ(3rgAw>9#Z@|I`vQttvpS$_V+2EiV+-o zBm^sknf%h};#kcbwniEpmUAngYK`R3Bq5j|Oq=+paiMRs*qz81_n)`rQewcN8X?$J zm^C-L-;LMLWR1&d#+s)mXrP2c$An-hFz2TwQ>V_uY}S>`Rm*fTseBZNT!df`j|g+- z-Zkwd!xGpivzaTKb5dx_Xby!6!CYZF?YNU2qqduAjL;byem{2nB6&- zOZ1-aChbYzxHZjPO-+nAR1>Ze#4lTzwck}3%zKqceq6Efy!C(t>OPi3JiiUC|s{`>vX( zJ(?4jD664GlR2a(1kr02W|`NZyxm?kL{d67zG`kARZQWKjS!>>=7*0aQ}V)V$;z%I zYm=OMsx_5E9zqcPPr{V+xppIKbsdr4)w%S`=>|%sap;r~#2V(tPp6IzZmT2PU;aGg zvC{<_Xv(4eLXh1sJrj@K+xe)D*hu{}SA<@m^64CM5Q6N5DSP;hbYm1R6#6en9Ab8X zwwQ6KR0xs@)4NJpdA+)pOf(wPrtof{N^=g47J_)f{PoYF_ZEfJ5cT83Y4M&7L(FnpI<1!a~TdM*TMg!!b&+;L-mA+e~-%yaiGp`BJ7S}6ok z!t`A*e#_H2a+2Fzbfl_JHkHigkY$TX5XXhW^cxcM+s8UT5@fk$zDb~zI?mxx=4_Q9 z`TD}Vf5>dVqcn{zH63KNd2$++5)QRaRS7l`=B#I9bk)5JS=A5LL1pQuX^{ekOzSwO`U3^OvEtmDkyqPVsEIv?8iio%VJ^UDE+F?M z9KtV{9n1t=(G<6fci@mKZkGss++lX%34Jt1UQYRH4voRfc^2kryvf(_nzyXs&=b7oH8AzDClq4C zP=3#$Dr^`dVSdKeQ-Q6=WF3bZvGsiXvM*tu>K6MJbzjdRB{s2Um|ocNzQcAH@6I7} zY==QG|Glv?4-S=MV@=1#8h~AR4Yue`PY$ia7A=Drg#A1eo4v$~Ll7hcrWsC#A{-lz z8#%;qY|Mi>7-x(>4j`#FhkS4V4TE_Vr_~G`Wko(5vc^%?4D&4xMPnR_3SSQC;!y03 zLs1*Yt1*sOtxX&njN`S(Qk=bC&SIR+l=yS#Fb?W7UruqH>eq4P2mZjJ^EmS5Fdv~l zcq-}zD&NeZMidXpFuPD`{Dm^4v0m@pe&>jT$b>wfJou>JftdiE=Gzk0|_t zri;RlZRgM+6n=MM`lCuzLkZcrgF~ZGLKeeZjC#`z#i(Q#hYC@Q=E1xhC(7DQC~F;e zb0`;OtpsK50@T3D5K$mYgE(Z50y!LJ5vu9eD7B0BaA^2TRo%@xV16mu``#8syJ9bg zrU}80qG;D-qWsT5`L7ksp){2LdN6OHSI9>rL4M>ALnCn;<}c_vlF^a`hH%IUElCc{ zJoGPLw}|FOF6Gc2G&eGs)6ogd6`c`n3FVLx+N9wy&!DI3qZAF7@+S^GK*J@2S&r_^ z60Mubehzh}tNJ+R0n^~MXacvP33LzRkQYbI>rrXC*u!r$O7$T156wAo&{)3 zSvZH#r80k*8t7u1(BgJRaHth6t}aX`^us66EK8y|)Qx624Q3rW>teLsj?o~7{gJ8O} zi+G|2@q|`9hwKqgWFejykSl@?TLTevkOU4nBk0gU(D4eXNM9onnFJo<&`?As0WfQj zuQ-Ujg~}7PdrMm*L@2WXW&uK*eF$;-A=Fugkf$A?j|SqcX9$IEBP7xc7V*|}gh;Lk zmE;JS79riL5Gj+JX^jBUk*GjXF7+_(xEY5@l3h=CR-i43(JVepp-H5M~Xz|8kz MC9<$GiCC%s14xJ&IRF3v diff --git a/python/tests/dispersive_eigenmode.py b/python/tests/dispersive_eigenmode.py index 60521b756..afe20b3f6 100644 --- a/python/tests/dispersive_eigenmode.py +++ b/python/tests/dispersive_eigenmode.py @@ -63,7 +63,7 @@ def verify_output_and_slice(self,material,frequency): sim.init_sim() # Check to make sure the get_slice routine is working with frequency - n_slice = np.sqrt(np.max(sim.get_epsilon(frequency))) + n_slice = np.sqrt(np.max(sim.get_epsilon(frequency=frequency,snap=True))) self.assertAlmostEqual(n,n_slice, places=4) # Check to make sure h5 output is working with frequency diff --git a/python/tests/faraday_rotation.py b/python/tests/faraday_rotation.py index 140868b5b..73485758c 100644 --- a/python/tests/faraday_rotation.py +++ b/python/tests/faraday_rotation.py @@ -41,12 +41,11 @@ def check_rotation(self, mat, L, fsrc, zsrc, resolution, tmax, zout, kpred, tol= boundary_layers=pml_layers, default_material=mat, resolution=resolution) - record_vol = mp.Volume(center=mp.Vector3(0, 0, zout)) record_Ex, record_Ey, record_t = [], [], [] def record_ex_ey(sim): - record_Ex.append(sim.get_array(vol=record_vol, component=mp.Ex)) - record_Ey.append(sim.get_array(vol=record_vol, component=mp.Ey)) + record_Ex.append(sim.get_field_point(mp.Ex, mp.Vector3(0, 0, zout))) + record_Ey.append(sim.get_field_point(mp.Ey, mp.Vector3(0, 0, zout))) record_t.append(sim.meep_time()) self.sim.run(mp.after_time(0.5*tmax, mp.at_every(1e-6, record_ex_ey)), until=tmax) diff --git a/python/tests/simulation.py b/python/tests/simulation.py index 883ce050c..b874b3eae 100644 --- a/python/tests/simulation.py +++ b/python/tests/simulation.py @@ -373,10 +373,10 @@ def test_get_array_output(self): mp.output_tot_pwr(sim) mp.output_efield(sim) - eps_arr = sim.get_epsilon() - efield_z_arr = sim.get_efield_z() - energy_arr = sim.get_tot_pwr() - efield_arr = sim.get_efield() + eps_arr = sim.get_epsilon(snap=True) + efield_z_arr = sim.get_efield_z(snap=True) + energy_arr = sim.get_tot_pwr(snap=True) + efield_arr = sim.get_efield(snap=True) fname_fmt = os.path.join(self.temp_dir, 'test_get_array_output-{}-000020.00.h5') diff --git a/src/array_slice.cpp b/src/array_slice.cpp index a156ba050..4d562ecde 100644 --- a/src/array_slice.cpp +++ b/src/array_slice.cpp @@ -76,6 +76,7 @@ typedef struct { component invmu_cs[3]; direction invmu_ds[3]; + bool snap; } array_slice_data; #define UNUSED(x) (void)x // silence compiler warnings @@ -332,7 +333,7 @@ static void get_array_slice_chunkloop(fields_chunk *fc, int ichnk, component cgr frequency)); if (abs(tr) == 0.0) tr += 4.0; // default inveps == 1 } - fields[i] = (4.0 * data->ninveps) / tr; + fields[i] = (data->snap ? 1.0 : IVEC_LOOP_WEIGHT(s0, s1, e0, e1, 1.0)) * (4.0 * data->ninveps) / tr; } else if (cS[i] == Permeability) { complex tr(0.0, 0.0); @@ -344,7 +345,7 @@ static void get_array_slice_chunkloop(fields_chunk *fc, int ichnk, component cgr frequency)); if (abs(tr) == 0.0) tr += 4.0; // default invmu == 1 } - fields[i] = (4.0 * data->ninvmu) / tr; + fields[i] = (data->snap ? 1.0 : IVEC_LOOP_WEIGHT(s0, s1, e0, e1, 1.0)) * (4.0 * data->ninvmu) / tr; } else { double f[2]; @@ -355,7 +356,7 @@ static void get_array_slice_chunkloop(fields_chunk *fc, int ichnk, component cgr fc->f[cS[i]][k][idx + off[2 * i] + off[2 * i + 1]]); else f[k] = 0; - fields[i] = complex(f[0], f[1]) * ph[i]; + fields[i] = (data->snap ? 1.0 : IVEC_LOOP_WEIGHT(s0, s1, e0, e1, 1.0)) * complex(f[0], f[1]) * ph[i]; } } @@ -470,44 +471,117 @@ int fields::get_array_slice_dimensions(const volume &where, size_t dims[3], dire return rank; } +/**********************************************************************/ +/* increment a multi-index in which n_i runs over*/ +/* the range 0 <= n_i < nMax_i. return true if this is the increment */ +/* that brings the multiindex to all zeros. */ +/**********************************************************************/ +bool increment(size_t *n, size_t *nMax, int rank) { + for (rank--; rank >= 0; rank--) + if (++n[rank] < nMax[rank]) + return false; + else + n[rank] = 0; + return true; +} + +// data_size = 1,2 for real,complex-valued array +double *collapse_array(double *array, int *rank, size_t dims[3], direction dirs[3], volume where, + int data_size = 1) { + + /*--------------------------------------------------------------*/ + /*- detect empty dimensions and compute rank and strides for */ + /*- collapsed array */ + /*--------------------------------------------------------------*/ + int full_rank = *rank; + if (full_rank == 0) return array; + + int reduced_rank = 0; + size_t reduced_dims[3], reduced_stride[3] = {1, 1, 1}; + direction reduced_dirs[3]; + for (int r = 0; r < full_rank; r++) { + if (where.in_direction(dirs[r]) == 0.0) + reduced_stride[r] = 0; // degenerate dimension, to be collapsed + else { + reduced_dirs[reduced_rank] = dirs[r]; + reduced_dims[reduced_rank++] = dims[r]; + } + } + if (reduced_rank==0) { + *rank = 0; + return array; // return array as is for singleton use case + } + if (reduced_rank == full_rank) return array; // nothing to collapse + + /*--------------------------------------------------------------*/ + /*- set up strides into full and reduced arrays */ + /*--------------------------------------------------------------*/ + size_t stride[3] = {1, 1, 1}; // non-reduced array strides + if (full_rank == 2) + stride[0] = dims[1]; // rstride is already all set in this case + else if (full_rank == 3) { + stride[0] = dims[1] * dims[2]; + stride[1] = dims[2]; + if (reduced_rank == 2) + reduced_stride[reduced_stride[0] != 0 ? 0 : 1] = reduced_dims[1]; + } + + /*--------------------------------------------------------------*/ + /*- allocate reduced array and compress full array into it -*/ + /*--------------------------------------------------------------*/ + size_t reduced_grid_size = reduced_dims[0] * (reduced_rank == 2 ? reduced_dims[1] : 1); + size_t reduced_array_size = data_size * reduced_grid_size; + double *reduced_array = new double[reduced_array_size]; + if (!reduced_array) abort("%s:%i: out of memory (%zu)", __FILE__, __LINE__, reduced_array_size); + memset(reduced_array, 0, reduced_array_size * sizeof(double)); + + size_t n[3] = {0, 0, 0}; + do { + size_t index = n[0] * stride[0] + n[1] * stride[1] + n[2] * stride[2]; + size_t rindex = n[0] * reduced_stride[0] + n[1] * reduced_stride[1] + n[2] * reduced_stride[2]; + for (int i = 0; i < data_size; i++) + reduced_array[data_size * rindex + i] += array[data_size * index + i]; + } while (!increment(n, dims, full_rank)); + + *rank = reduced_rank; + for (int r = 0; r < reduced_rank; r++) { + dims[r] = reduced_dims[r]; + dirs[r] = reduced_dirs[r]; + } + delete[] array; + return reduced_array; +} + +cdouble *collapse_array(cdouble *array, int *rank, size_t dims[3], direction dirs[3], + volume where) { + return (cdouble *)collapse_array((double *)array, rank, dims, dirs, where, 2); +} + /**********************************************************************/ /* precisely one of fun, rfun, should be non-NULL */ /**********************************************************************/ void *fields::do_get_array_slice(const volume &where, std::vector components, field_function fun, field_rfunction rfun, void *fun_data, - void *vslice, double frequency) { + void *vslice, double frequency, bool snap) { am_now_working_on(FieldOutput); /***************************************************************/ /* call get_array_slice_dimensions to get slice dimensions and */ /* partially initialze an array_slice_data struct */ /***************************************************************/ - // by tradition, empty dimensions in time-domain field arrays are *not* collapsed - // TODO make this a caller-specifiable parameter to get_array_slice()? - bool collapse = false, snap = true; size_t dims[3]; direction dirs[3]; array_slice_data data; - get_array_slice_dimensions(where, dims, dirs, collapse, snap, 0, &data); + int rank = get_array_slice_dimensions(where, dims, dirs, false, snap, 0, &data); size_t slice_size = data.slice_size; - bool complex_data = (rfun == 0); - cdouble *zslice; - double *slice; - if (vslice == 0) { - if (complex_data) { - zslice = new cdouble[slice_size]; - memset(zslice, 0, slice_size * sizeof(cdouble)); - vslice = (void *)zslice; - } - else { - slice = new double[slice_size]; - memset(slice, 0, slice_size * sizeof(double)); - vslice = (void *)slice; - } - } + int elem_size = complex_data ? 2 : 1; + void *vslice_uncollapsed; - data.vslice = vslice; + vslice_uncollapsed = memset(new double[slice_size * elem_size], 0, slice_size * elem_size * sizeof(double)); + + data.vslice = vslice_uncollapsed; + data.snap = snap; data.fun = fun; data.rfun = rfun; data.fun_data = fun_data; @@ -550,12 +624,22 @@ void *fields::do_get_array_slice(const volume &where, std::vector com ++data.ninvmu; } - loop_in_chunks(get_array_slice_chunkloop, (void *)&data, where, Centered, true, true); + loop_in_chunks(get_array_slice_chunkloop, (void *)&data, where, Centered, true, snap); - /***************************************************************/ - /*consolidate full array on all cores */ - /***************************************************************/ - vslice = (void *)array_to_all((double *)vslice, (complex_data ? 2 : 1) * slice_size); + if (!snap) { + double *slice = collapse_array((double *)vslice_uncollapsed, &rank, dims, dirs, where, elem_size); + rank = get_array_slice_dimensions(where, dims, dirs, true, false, 0, &data); + slice_size = data.slice_size; + vslice_uncollapsed = (double*) slice; + } + if (vslice) { + memcpy(vslice, vslice_uncollapsed, slice_size * elem_size * sizeof(double)); + delete[] (double*) vslice_uncollapsed; + } + else + vslice = vslice_uncollapsed; + + array_to_all((double *)vslice, elem_size * slice_size); delete[] data.offsets; delete[] data.fields; @@ -571,168 +655,93 @@ void *fields::do_get_array_slice(const volume &where, std::vector com /***************************************************************/ double *fields::get_array_slice(const volume &where, std::vector components, field_rfunction rfun, void *fun_data, double *slice, - double frequency) { + double frequency, bool snap) { return (double *)do_get_array_slice(where, components, 0, rfun, fun_data, (void *)slice, - frequency); + frequency, snap); } cdouble *fields::get_complex_array_slice(const volume &where, std::vector components, field_function fun, void *fun_data, cdouble *slice, - double frequency) { + double frequency, bool snap) { return (cdouble *)do_get_array_slice(where, components, fun, 0, fun_data, (void *)slice, - frequency); + frequency, snap); } -double *fields::get_array_slice(const volume &where, component c, double *slice, double frequency) { +double *fields::get_array_slice(const volume &where, component c, double *slice, double frequency, + bool snap) { std::vector components(1); components[0] = c; return (double *)do_get_array_slice(where, components, 0, default_field_rfunc, 0, (void *)slice, - frequency); + frequency, snap); } double *fields::get_array_slice(const volume &where, derived_component c, double *slice, - double frequency) { + double frequency, bool snap) { int nfields; component carray[12]; field_rfunction rfun = derived_component_func(c, gv, nfields, carray); std::vector cs(carray, carray + nfields); - return (double *)do_get_array_slice(where, cs, 0, rfun, &nfields, (void *)slice, frequency); + return (double *)do_get_array_slice(where, cs, 0, rfun, &nfields, (void *)slice, + frequency, snap); } cdouble *fields::get_complex_array_slice(const volume &where, component c, cdouble *slice, - double frequency) { + double frequency, bool snap) { std::vector components(1); components[0] = c; return (cdouble *)do_get_array_slice(where, components, default_field_func, 0, 0, (void *)slice, - frequency); + frequency, snap); } cdouble *fields::get_source_slice(const volume &where, component source_slice_component, cdouble *slice) { - size_t dims[3]; direction dirs[3]; vec min_max_loc[2]; - bool collapse = false, snap = false; - int rank = get_array_slice_dimensions(where, dims, dirs, collapse, snap, min_max_loc); + int rank = get_array_slice_dimensions(where, dims, dirs, false, false, min_max_loc); size_t slice_size = dims[0] * (rank >= 2 ? dims[1] : 1) * (rank == 3 ? dims[2] : 1); source_slice_data data; data.source_component = source_slice_component; data.slice_imin = gv.round_vec(min_max_loc[0]); data.slice_imax = gv.round_vec(min_max_loc[1]); - data.slice = slice ? slice : new cdouble[slice_size]; + data.slice = new cdouble[slice_size]; if (!data.slice) abort("%s:%i: out of memory (%zu)", __FILE__, __LINE__, slice_size); - loop_in_chunks(get_source_slice_chunkloop, (void *)&data, where, Centered, true, true); - return array_to_all(data.slice, slice_size); -} - -/**********************************************************************/ -/* increment a multi-index in which n_i runs over*/ -/* the range 0 <= n_i < nMax_i. return true if this is the increment */ -/* that brings the multiindex to all zeros. */ -/**********************************************************************/ -bool increment(size_t *n, size_t *nMax, int rank) { - for (rank--; rank >= 0; rank--) - if (++n[rank] < nMax[rank]) - return false; - else - n[rank] = 0; - return true; -} + loop_in_chunks(get_source_slice_chunkloop, (void *)&data, where, Centered, true, false); -// data_size = 1,2 for real,complex-valued array -double *collapse_array(double *array, int *rank, size_t dims[3], direction dirs[3], volume where, - int data_size = 1) { - - /*--------------------------------------------------------------*/ - /*- detect empty dimensions and compute rank and strides for */ - /*- collapsed array */ - /*--------------------------------------------------------------*/ - int full_rank = *rank; - if (full_rank == 0) return array; - - int reduced_rank = 0; - size_t reduced_dims[3], reduced_stride[3] = {1, 1, 1}; - direction reduced_dirs[3]; - for (int r = 0; r < full_rank; r++) { - if (where.in_direction(dirs[r]) == 0.0) - reduced_stride[r] = 0; // degenerate dimension, to be collapsed - else { - reduced_dirs[reduced_rank] = dirs[r]; - reduced_dims[reduced_rank++] = dims[r]; - } - } - if (reduced_rank==0) { - *rank = 0; - return array; // return array as is for singleton use case - } - if (reduced_rank == full_rank) return array; // nothing to collapse + cdouble *slice_collapsed = collapse_array(data.slice, &rank, dims, dirs, where); + rank = get_array_slice_dimensions(where, dims, dirs, true, false); + slice_size = dims[0] * (rank >= 2 ? dims[1] : 1) * (rank == 3 ? dims[2] : 1); - /*--------------------------------------------------------------*/ - /*- set up strides into full and reduced arrays */ - /*--------------------------------------------------------------*/ - size_t stride[3] = {1, 1, 1}; // non-reduced array strides - if (full_rank == 2) - stride[0] = dims[1]; // rstride is already all set in this case - else if (full_rank == 3) { - stride[0] = dims[1] * dims[2]; - stride[1] = dims[2]; - if (reduced_rank == 2) - reduced_stride[reduced_stride[0] != 0 ? 0 : 1] = reduced_dims[1]; + if (slice != 0) { + for (size_t i = 0; i < slice_size; ++i) + slice[i] = slice_collapsed[i]; + slice = array_to_all(slice, slice_size); } - /*--------------------------------------------------------------*/ - /*- allocate reduced array and compress full array into it -*/ - /*--------------------------------------------------------------*/ - size_t reduced_grid_size = reduced_dims[0] * (reduced_rank == 2 ? reduced_dims[1] : 1); - size_t reduced_array_size = data_size * reduced_grid_size; - double *reduced_array = new double[reduced_array_size]; - if (!reduced_array) abort("%s:%i: out of memory (%zu)", __FILE__, __LINE__, reduced_array_size); - memset(reduced_array, 0, reduced_array_size * sizeof(double)); - - size_t n[3] = {0, 0, 0}; - do { - size_t index = n[0] * stride[0] + n[1] * stride[1] + n[2] * stride[2]; - size_t rindex = n[0] * reduced_stride[0] + n[1] * reduced_stride[1] + n[2] * reduced_stride[2]; - for (int i = 0; i < data_size; i++) - reduced_array[data_size * rindex + i] += array[data_size * index + i]; - } while (!increment(n, dims, full_rank)); - - *rank = reduced_rank; - for (int r = 0; r < reduced_rank; r++) { - dims[r] = reduced_dims[r]; - dirs[r] = reduced_dirs[r]; - } - delete[] array; - return reduced_array; -} + slice_collapsed = array_to_all(slice_collapsed, slice_size); -cdouble *collapse_array(cdouble *array, int *rank, size_t dims[3], direction dirs[3], - volume where) { - return (cdouble *)collapse_array((double *)array, rank, dims, dirs, where, 2); + return slice_collapsed; } /***************************************************************/ /***************************************************************/ /***************************************************************/ -std::vector fields::get_array_metadata(const volume &where, bool collapse_empty_dimensions, - bool snap_empty_dimensions) { +std::vector fields::get_array_metadata(const volume &where) { /* get extremal corners of subgrid and array of weights, collapsed if necessary */ size_t dims[3]; direction dirs[3]; vec min_max_loc[2]; // extremal points in subgrid - int rank = get_array_slice_dimensions(where, dims, dirs, false /*collapse_empty_dimensions*/, - snap_empty_dimensions, min_max_loc); + int rank = get_array_slice_dimensions(where, dims, dirs, true, false, min_max_loc); + int full_rank = rank; direction full_dirs[3]; for (int fr = 0; fr < rank; fr++) full_dirs[fr] = dirs[fr]; double *weights = get_array_slice(where, NO_COMPONENT); - if (collapse_empty_dimensions) weights = collapse_array(weights, &rank, dims, dirs, where); /* get length and endpoints of x,y,z tics arrays */ size_t nxyz[3] = {1, 1, 1}; @@ -740,7 +749,7 @@ std::vector fields::get_array_metadata(const volume &where, bool collaps for (int fr = 0, rr = 0; fr < full_rank; fr++) { direction d = full_dirs[fr]; int nd = d - X; - if (where.in_direction(d) == 0.0 && collapse_empty_dimensions) { + if (where.in_direction(d) == 0.0) { xyzmin[nd] = xyzmax[nd] = where.in_direction_min(d); nxyz[nd] = 1; } diff --git a/src/meep.hpp b/src/meep.hpp index 2a489f2ac..716f3a513 100644 --- a/src/meep.hpp +++ b/src/meep.hpp @@ -1640,13 +1640,8 @@ class fields { // in the x direction (corresponding to the two grid points // nearest x0, from which fields at x0 are interpolated). // if collapse_empty_dimensions==true, all such length-2 - // array dimensions are collaped to length 1 by doing the + // array dimensions are collapsed to length 1 by doing the // interpolation before returning the array. - // currently, collapse_empty_dimensions is always false for the - // time-domain arrays returned by get_field_array and always - // true for the frequency-domain arrays returned by get_dft_array, - // so an alternative name for `collapse_empty_dimensions` would be - // `is_dft_array`. // // the `data` parameter is used internally in get_array_slice // and should be ignored by external callers. @@ -1655,10 +1650,6 @@ class fields { bool snap_empty_dimensions = false, vec *min_max_loc = NULL, void *data = 0, component cgrid = Centered); - int get_dft_array_dimensions(const volume &where, size_t dims[3], direction dirs[3]) { - return get_array_slice_dimensions(where, dims, dirs, true); - } - // given a subvolume, return a column-major array containing // the given function of the field components in that subvolume // if slice is non-null, it must be a user-allocated buffer @@ -1667,24 +1658,27 @@ class fields { // must eventually be caller-deallocated via delete[]. double *get_array_slice(const volume &where, std::vector components, field_rfunction rfun, void *fun_data, double *slice = 0, - double frequency = 0); + double frequency = 0, + bool snap = false); std::complex *get_complex_array_slice(const volume &where, std::vector components, field_function fun, void *fun_data, std::complex *slice = 0, - double frequency = 0); + double frequency = 0, + bool snap = false); // alternative entry points for when you have no field // function, i.e. you want just a single component or // derived component.) double *get_array_slice(const volume &where, component c, double *slice = 0, - double frequency = 0); + double frequency = 0, bool snap = false); double *get_array_slice(const volume &where, derived_component c, double *slice = 0, - double frequency = 0); + double frequency = 0, bool snap = false); std::complex *get_complex_array_slice(const volume &where, component c, std::complex *slice = 0, - double frequency = 0); + double frequency = 0, + bool snap = false); // like get_array_slice, but for *sources* instead of fields std::complex *get_source_slice(const volume &where, component source_slice_component, @@ -1693,13 +1687,11 @@ class fields { // master routine for all above entry points void *do_get_array_slice(const volume &where, std::vector components, field_function fun, field_rfunction rfun, void *fun_data, void *vslice, - double frequency = 0); + double frequency = 0, bool snap = false); - /* fetch and return coordinates and integration weights of grid points covered by an array slice, - */ + /* fetch and return coordinates and integration weights of grid points covered by an array slice, */ /* packed into a vector with format [NX, xtics[:], NY, ytics[:], NZ, ztics[:], weights[:] ] */ - std::vector get_array_metadata(const volume &where, bool collapse_empty_dimensions = true, - bool snap_empty_dimensions = false); + std::vector get_array_metadata(const volume &where); // step.cpp methods: double last_step_output_wall_time; diff --git a/src/meepgeom.cpp b/src/meepgeom.cpp index 82e910353..c4860a009 100644 --- a/src/meepgeom.cpp +++ b/src/meepgeom.cpp @@ -2468,12 +2468,11 @@ void material_grids_addgradient(meep::realnum *v, size_t ng, std::complexget_array_slice_dimensions(where, &dims[3 * c], dirs, collapse, snap, min_max_loc, + int rank = f->get_array_slice_dimensions(where, &dims[3 * c], dirs, true, false, min_max_loc, 0, field_dir[c]); vector3 max_corner = vec_to_vector3(min_max_loc[1]); diff --git a/tests/array-metadata.cpp b/tests/array-metadata.cpp index 3c22b0c18..1967ae719 100644 --- a/tests/array-metadata.cpp +++ b/tests/array-metadata.cpp @@ -68,15 +68,15 @@ bool equal_float(double *array1, double *array2, int N) { /* if the environment variable MEEP_ARRAY_METADATA_LOGFILE is */ /* set, more detailed output is written to that file. */ /***************************************************************/ -bool test_array_metadata(meep::fields &f, const volume &where, bool collapse_empty_dimensions) { +bool test_array_metadata(meep::fields &f, const volume &where) { /***************************************************************/ /* step 1: get coordinate grids and weights as reported by */ /* get_array_metadata */ /***************************************************************/ size_t dims[3]; direction dirs[3]; - int rank = f.get_array_slice_dimensions(where, dims, dirs, collapse_empty_dimensions); - std::vector xyzw = f.get_array_metadata(where, collapse_empty_dimensions); + int rank = f.get_array_slice_dimensions(where, dims, dirs); + std::vector xyzw = f.get_array_metadata(where); // convert to a more convenient format int offset = 0; @@ -225,8 +225,6 @@ int main(int argc, char *argv[]) { double res = 10.0; - bool collapse_empty_dimensions = false; - // double-valued command-line parameters vector parm_name; vector parm_adrs; @@ -255,13 +253,6 @@ int main(int argc, char *argv[]) { /*- parse arguments --------------------------------------------*/ /*--------------------------------------------------------------*/ for (int narg = 1; narg < argc; narg++) { - // process boolean-valued parameters - if (!strcasecmp(argv[narg], "--collapse")) { - collapse_empty_dimensions = true; - printf("Collapsing empty array dimensions.\n"); - continue; - } - // process double-valued parameters size_t np; for (np = 0; np < parm_name.size(); np++) @@ -327,6 +318,6 @@ int main(int argc, char *argv[]) { vmax.set_direction(Z, vzmax); } volume slice(vmin, vmax); - bool test_passed = test_array_metadata(f, slice, collapse_empty_dimensions); + bool test_passed = test_array_metadata(f, slice); return test_passed ? 0 : 1; } diff --git a/tests/array-slice-ll-ref.h5 b/tests/array-slice-ll-ref.h5 index 9d359db83c4e00f0f26c222cc99daff40af22f63..c54e9c945a0438157e822e0d91ce0f91332f5084 100644 GIT binary patch delta 53 zcmaEGn(4u5rVVSDS*9|-iQl||IfMz!n!JIzVX^|V(BxTcUF;w!28N7{8%-HECvY^( F1pra=6tMsR delta 53 zcmaEGn(4u5rVVSDSvF2S8@YJ{a|jcdHF*Pb!(;_!p~