inconsistent DFT fields monitors of the adjoint solver

[oskooi]

Currently, the EigenmodeCoefficient class of the adjoint solver defines the DFT field monitors on the Yee grid:

meep/python/adjoint/objective.py

Lines 154 to 162 in 434f3fb

	def register_monitors(self, frequencies):
	self._frequencies = np.asarray(frequencies)
	self._monitor = self.sim.add_mode_monitor(
	frequencies,
	mp.ModeRegion(center=self.volume.center, size=self.volume.size),
	yee_grid=True,
	decimation_factor=self.decimation_factor,
	)
	return self._monitor

However, the FourierFields class defines the DFT field monitors at the voxel centers:

meep/python/adjoint/objective.py

Lines 239 to 248 in 434f3fb

	def register_monitors(self, frequencies):
	self._frequencies = np.asarray(frequencies)
	self._monitor = self.sim.add_dft_fields(
	[self.component],
	self._frequencies,
	where=self.volume,
	yee_grid=False,
	decimation_factor=self.decimation_factor,
	)
	return self._monitor

And yet, the Near2FarFields class defines the DFT field monitors on the Yee grid:

meep/python/adjoint/objective.py

Lines 340 to 347 in 434f3fb

	def register_monitors(self, frequencies):
	self._frequencies = np.asarray(frequencies)
	self._monitor = self.sim.add_near2far(
	self._frequencies,
	*self.Near2FarRegions,
	decimation_factor=self.decimation_factor,
	)
	return self._monitor

Note: for Near2FarFields, the yee_grid option is not defined here but use_centered_grid (third last argument) is set to false in fields::add_near2far:

meep/src/near2far.cpp

Line 635 in 849bf3d

F = add_dft(c, w->v, freq, Nfreq, true, s * w->weight, F, false, 1.0, false, c0, decimation_factor);

meep/src/meep.hpp

Lines 1961 to 1966 in 434f3fb

	dft_chunk *add_dft(component c, const volume &where, double freq_min, double freq_max, int Nfreq,
	bool include_dV_and_interp_weights = true,
	std::complex<double> stored_weight = 1.0, dft_chunk *chunk_next = 0,
	bool sqrt_dV_and_interp_weights = false,
	std::complex<double> extra_weight = 1.0, bool use_centered_grid = true,
	int vc = 0, int decimation_factor = 0) {

When defining an objective function as part of an adjoint simulation, the user therefore needs to be consistent with these definitions. This is similar to the set up used in the unit test for the adjoint solver where the yee_grid option is explicitly defined:

meep/python/tests/test_adjoint_solver.py

Lines 79 to 91 in 434f3fb

	if mon_type.name == 'EIGENMODE':
	mode = sim.add_mode_monitor(frequencies,
	mp.ModeRegion(center=mp.Vector3(0.5*sxy-dpml-0.1),
	size=mp.Vector3(0,sxy-2*dpml,0)),
	yee_grid=True,
	eig_parity=eig_parity)
	elif mon_type.name == 'DFT':
	mode = sim.add_dft_fields([mp.Ez],
	frequencies,
	center=mp.Vector3(1.25),
	size=mp.Vector3(0.25,1,0),
	yee_grid=False)

Ensuring this type of consistency (particularly for the mode coefficients which does not use the default option) can be problematic and a source of error for the average user. It would be useful if the DFT fields monitors for all objective classes were always defined at the voxel centers.

[ahoenselaar]

Beyond consistency, are there reasons to favor monitors on the centered grid over monitors on the Yee grid?

The current implementation of monitors on the centered grid is inefficient because it performs interpolation whenever the monitor data is updated during time stepping. It would be more efficient to place DFT monitors on the Yee grid and conditionally perform interpolation when returning data to the user. This is supported by a linearity argument for the DFT and the interpolation step.

[stevengj]

The basic idea here is that the centered grid is much easier for the user to deal with than the Yee grid. And FourierFields are designed to be used directly by the user. e.g. it would be pretty hard for the user to compute ExH if all the field components are at different points in space.

However, it would be good to have the option to get FourierFields directly on the Yee grid, e.g. the version that @mawc2019 is working on.

Whereas the DFT fields for the eigenmode coefficients and the near2far transformations are mainly for use by Meep, and Meep knows how to handle Yee-grid data.

It would be more efficient to place DFT monitors on the Yee grid and conditionally perform interpolation when returning data to the user.

I agree: it would be mathematically equivalent to do the DFT accumulation on the Yee grid and only convert to the centered grid when we output to the user.

[oskooi]

Setting yee_grid=False for the EigenModeCoefficient class of objective.py is causing two of the gradient tests in test_adjoint_solver.py to fail in a significant way:

======================================================================
FAIL: test_adjoint_solver_eigenmode (__main__.TestAdjointSolver)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "./tests/test_adjoint_solver.py", line 320, in test_adjoint_solver_eigenmode
    self.assertClose(adjsol_obj,S12_unperturbed,epsilon=1e-6)
  File "/home/install/meep/python/tests/utils.py", line 28, in assertClose
    self.assertLessEqual(diff_norm, epsilon * np.maximum(x_norm, y_norm), msg)
AssertionError: 0.8343982855925702 not less than or equal to 3.0114380814839494e-05 : 

======================================================================
FAIL: test_gradient_backpropagation (__main__.TestAdjointSolver)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "./tests/test_adjoint_solver.py", line 362, in test_gradient_backpropagation
    self.assertClose(adjsol_obj,S12_unperturbed,epsilon=1e-6)
  File "/home/install/meep/python/tests/utils.py", line 28, in assertClose
    self.assertLessEqual(diff_norm, epsilon * np.maximum(x_norm, y_norm), msg)
AssertionError: 0.8266286252600565 not less than or equal to 2.9548554303440535e-05 : 

----------------------------------------------------------------------
Ran 4 tests in 131.852s

FAILED (failures=2)

@smartalecH, why was it necessary to override the yee_grid=False default (see below) of add_mode_monitor for EigenModeCoefficient? Is it just a matter of adding the interpolation step from the Yee to the centered grid to the backpropagation step when computing the gradients to enable using the centered grid for EigenModeCoefficient?

meep/python/simulation.py

Line 3027 in c65a149

yee_grid = kwargs.get("yee_grid", False)

[smartalecH]

Yes that's expected... we must evaluate on the yee grid.

[smartalecH]

(...to continue using a flipped eigenmode source as our adjoint source. We can restrict the fields, in theory, if we store the field profiles at the monitors and process those accordingly, just like we do with the Near2Far.

Although we still need to finish #1547 before we can proceed with that approach.)