User-defined sources and mediums

[tylerflex]

[momchil-flex]

We need to put a bit more thought in how to set this up, keeping in mind the physics of what should happen. First of all, maybe we want to do just CustomFieldSource for now rather than also custom current source, which is less likely to be needed. One point related to that is that this seems a bit fishy to me. I feel like instead of doing that we should just have a CurrentData or something.

Regarding CustomFieldSource - I think this only makes sense as a 2D object. Basically you use equivalence principle to inject J from the provided H and M from the provided E. I also checked how importing sources works in Lumerical, and I have a few things I'm wondering about:

• In Lumerical, you don't specify source geometry - the position along the normal direction, and the size of the source plane, is taken from the supplied data. Do we want to do something like that, or do we still require center and size, require that the source is planar, and validate that the supplied data covers the source geometry? But more precisely we should first recenter the data to the source center I would think, such that I can provide the fields computed doesn't matter where, and then they'll be inserted in the right source location (is this a frontend or backend thing to do?)
• In Lumerical, the E fields are mandatory, and the H fields are optional. If they are not provided, they do a numerical computation of the H fields from the E fields, which is only correct for smoothly varying modes - but they do that to attempt to still have the source being directional. If we instead use only the E fields to inject currents, the source will be injected in both directions, and normalization won't work. This is more of a backend implementation issue I guess but we need to decide how we want to handle it.
• direction should still be a Field of the source, so that it can be injected backward if needed.

[tylerflex]

... maybe we want to do just CustomFieldSource for now rather than also custom current source, which is less likely to be needed. One point related to that is that this seems a bit fishy to me. I feel like instead of doing that we should just have a CurrentData or something.

I agree with this, actually. We could just remove CustomCurrentSource for now and make the source accept a FieldData with .monitor set to None (changing MonitorData to make this optional?)

Regarding CustomFieldSource - I think this only makes sense as a 2D object. Basically you use equivalence principle to inject J from the provided H and M from the provided E.

Agree with this, we can assert_plane for this source.

In Lumerical, you don't specify source geometry - the position along the normal direction, and the size of the source plane, is taken from the supplied data. Do we want to do something like that, or do we still require center and size, require that the source is planar, and validate that the supplied data covers the source geometry?

Hm, I think the user specifying center and size could be nice in that it might tell us precisely where the user wants to inject the source. For example, we could have a rule similar to the monitor, but in reverse: user supplies the data and some volume specified by the source geometry. The extent of the data must be sufficiently large that we can interpolate the J, M into all the yee cells within the volume?

But more precisely we should first recenter the data to the source center I would think, such that I can provide the fields computed doesn't matter where, and then they'll be inserted in the right source location (is this a frontend or backend thing to do?)

I dont know if I agree with this. So you're saying that a (x,y,z) == (0,0,0) in the field_data would correspond to that being injected at the monitor center? I imagine this could be a headache. For example, let's say I measured fields using a FieldMonitor and then want to create a custom source using the FieldData (field_data). I might do something (using the current API) like

center = field_data.monitor.center
size = field_data.monitor.size

src = CustomFieldSource(center=center, size=size, direction="+", **field_data.field_components)

and not have to worry about reentering the field_data to the origin.

In Lumerical, the E fields are mandatory, and the H fields are optional. If they are not provided, they do a numerical computation of the H fields from the E fields, which is only correct for smoothly varying modes - but they do that to attempt to still have the source being directional. If we instead use only the E fields to inject currents, the source will be injected in both directions, and normalization won't work. This is more of a backend implementation issue I guess but we need to decide how we want to handle it. Direction should still be a Field of the source, so that it can be injected backward if needed.

Hm, seems like one issue is whether we want to split this into a directional source or just a literal source. Perhaps how we could do it is as follows:

• Add an optional direction field to the CurrentSource.
• If direction is unset, just inject literally what is supplied, with no attempt to make it directional.
• If direction is set, process the supplied field components (probably on the backend) to generate a directional source. (remove normal components, set tangential H as necessary).

Does this make some sense?

To summarize:

1. How do we specify the data? Whether to keep CurrentSource.Ex etc. as ScalarFieldDataArray fields or add a FieldData field and modify FieldData to generalize it to cases where we mightn't need a monitor. Advantage of current implementation is that we keep monitors in MonitorData, but it also means we reproduce code and is limiting in cases.
2. Is a (x,y,z) in the data relative to the monitor center? or is it in global coordinates? I like global coordinates because the FieldData returned by solver is already in global coordinates, making conversion easy.
3. Do we allow directional and non directional source?

[momchil-flex]

2. Is a (x,y,z) in the data relative to the monitor center? or is it in global coordinates? I like global coordinates because the FieldData returned by solver is already in global coordinates, making conversion easy.

I don't think I follow your reasoning, or maybe mine wasn't very clear. The problem with global coordinates is that it is only easy when you want to place a source exactly where a monitor used to be. Say that I have some simulation and a yz monitor at x = 10. Then I want to take those fields and inject them in another simulation but at x = -5. If we expect the FieldData provided to the CustomSource to cover the source geometry (center, size) region, then the user has to either shift all of the FieldData x coordinates manually, or shift their whole simulation so that the source is at x = 10. What I'm suggesting instead is that the FieldData is taken such that the center along x, y, z, is placed at the source center. That doesn't mean the (0, 0, 0) of the FieldData is placed at the monitor center, but rather we can compute the bounds of the provided data based on the present coordinates, and we can define the data center based on those bounds. We put that center at the source center. So it will directly work for the use case I just describe (as well as when the source is exactly where the monitor used to be). Similarly any offset in yz would just work right off the bat.

3. Do we allow directional and non directional source?

If direction is unset, just inject literally what is supplied, with no attempt to make it directional. If direction is set, process the supplied field components (probably on the backend) to generate a directional source. (remove normal components, set tangential H as necessary).

I don't understand how this would work. We inject currents, not fields. The only way I can see it happening is through the equivalence principle for the surface of the source, so normal components have to be removed.

Then, if only E is supplied and we do nothing, this would inject in both directions. Which could be fine, the question here is whether we want to replicate more closely the Lumerical behavior or not.

[tylerflex]

I see your point about the centering of the monitor. I think I missed the bit about the center of the data being placed at the source center. I instead thought you meant that the data coordinates were relative to the source center. I think what you said makes sense and the use case for placing a monitor at a different position (ie monitor at x=10, source at x=-5) is reasonable.

I don't understand how this would work. We inject currents, not fields. The only way I can see it happening is through the equivalence principle for the surface of the source, so normal components have to be removed.

I was under the impression we could just substitute J and M currents in a volume for a given E, H field, similar to how we convert a field monitor into an adjoint source. But we can limit the scope to have the user specify tangential E (or H) and we'll compute th J, M surface currents using the equivalence principle.

Then, if only E is supplied and we do nothing, this would inject in both directions. Which could be fine, the question here is whether we want to replicate more closely the Lumerical behavior or not.

Can you link the lumerical docs for this? I think if users will always want to use this to inject equivalent sources, we can just make direction mandatory and do some processing to make it unidirectional. However, I wonder if there are applications where the "raw" (propagating in both directions) fields are desired.

Ok to summarize again where we are here: when figuring out the J(r) to inject given a CustomFIeldSource with some field_data info E(r), on the front end:

• Make sure the data bounds are all large enough to interpolate into the volume defined by the source.bounds. on the back end:
• shift E(r) in the data to E(r-r') where r' is the "center" of the supplied data, as defined via the average of the bounds in x,y,z.
• Remove any normal field components, defined by the dimension where source.size == 0.
• Compute J(r) and M(r) from tangential components of E(r) and H(r).

There's still an open question of how to handle the frequency-dependence. We have source.source_time and field_data.f. A few options:

1. Ignore field_data.f or make a new scalar field data array that omits it. Use the source.source_time to define the J(t) and the field_data to define the J(r).
2. Allow only one field_data.f and enforce that it matches source_time.freq0 or at least fits within the frequency_range of the source.
3. Use field_data.f to define the source_time.freq0. If many freqs present, perhaps generate a superposition of current sources each with a different source_time.freq0. (How to chose fwidth?) where one might want to inject several J(r) patterns at different frequencies.

[momchil-flex]

I see your point about the centering of the monitor. I think I missed the bit about the center of the data being placed at the source center. I instead thought you meant that the data coordinates were relative to the source center. I think what you said makes sense and the use case for placing a monitor at a different position (ie monitor at x=10, source at x=-5) is reasonable.

I don't understand how this would work. We inject currents, not fields. The only way I can see it happening is through the equivalence principle for the surface of the source, so normal components have to be removed.

I was under the impression we could just substitute J and M currents in a volume for a given E, H field, similar to how we convert a field monitor into an adjoint source. But we can limit the scope to have the user specify tangential E (or H) and we'll compute th J, M surface currents using the equivalence principle.

Then, if only E is supplied and we do nothing, this would inject in both directions. Which could be fine, the question here is whether we want to replicate more closely the Lumerical behavior or not.

Can you link the lumerical docs for this? I think if users will always want to use this to inject equivalent sources, we can just make direction mandatory and do some processing to make it unidirectional. However, I wonder if there are applications where the "raw" (propagating in both directions) fields are desired.

The way this last aspect is handled in Lumerical is that if the user wants to explicitly use the E-field only, they just provide H fields with all values set to zero (as opposed to not providing an H field). https://optics.ansys.com/hc/en-us/articles/360034383014-Import-source-Simulation-object

Ok to summarize again where we are here: when figuring out the J(r) to inject given a CustomFIeldSource with some field_data info E(r), on the front end:

• Make sure the data bounds are all large enough to interpolate into the volume defined by the source.bounds. on the back end:

I think the ones below could be restricted to the backend?

• shift E(r) in the data to E(r-r') where r' is the "center" of the supplied data, as defined via the average of the bounds in x,y,z.
• Remove any normal field components, defined by the dimension where source.size == 0.
• Compute J(r) and M(r) from tangential components of E(r) and H(r).

===

There's still an open question of how to handle the frequency-dependence. We have source.source_time and field_data.f. A few options:

1. Ignore field_data.f or make a new scalar field data array that omits it. Use the source.source_time to define the J(t) and the field_data to define the J(r).
2. Allow only one field_data.f and enforce that it matches source_time.freq0 or at least fits within the frequency_range of the source.
3. Use field_data.f to define the source_time.freq0. If many freqs present, perhaps generate a superposition of current sources each with a different source_time.freq0. (How to chose fwidth?) where one might want to inject several J(r) patterns at different frequencies.

I think what we should support right now is fields at a single frequency. Your suggestion 3. is a general broadband source improvement that should be added to all sources, not specifically to the CustomFieldSource case. I think in the current setup it is mandatory for the user to set a frequency to the ScalarFieldData? This is not ideal as in the simplest case people just want to provide the spatial data, but it's maybe not too much to ask them to put the frequency in too. I think what we should simply enforce is that there's only one f but not that it matches source_time.freq0. Instead it's assumed the fields are "frequency independent". If there's more than one f a NotImplementedError for multi-frequency sources should be raised.

[tylerflex]

Sounds good to me overall

=========

I forgot to add a note that to the processing (centering, normal field removal, etc) should be done on backend, what you wrote makes sense there.

Only thing we should do front end is validate:

• The data bounds are all outside of the source geometry bounds.
• Make sure at least 1 tangential field is supplied.
• Make sure if we are doing uni-directional field (based on inputs) that the "direction" is supplied, more on this below.

=========

Ok so allow only one scalar_field.f value and maybe just warn if it's outside of the source.source_time.frequency_range? Sounds reasonable to me.

We can leave broadband to later, agreed.

=========

Regarding directional source:

So as far as I understand, if E and H are supplied, numerical will inject the equivalent surface currents "as is" (without processing them to make them uni-directional). But if H are omitted, the H are decided to make it unidirectional?

I think this could work. More specifically, I guess we would make it directional if both tangential H fields are omitted and require the direction to be set.

Do we want to handle if H is supplied but E is not?

[momchil-flex]

Ok so allow only one scalar_field.f value and maybe just warn if it's outside of the source.source_time.frequency_range? Sounds reasonable to me.

Yeah that sounds fine.

=========

Regarding directional source:

So as far as I understand, if E and H are supplied, numerical will inject the equivalent surface currents "as is" (without processing them to make them uni-directional). But if H are omitted, the H are decided to make it unidirectional?

I think this could work. More specifically, I guess we would make it directional if both tangential H fields are omitted and require the direction to be set.

Do we want to handle if H is supplied but E is not?

I think the direction parameter may have an effect even if the H fields are specified. It would just be adding a - sign to them or something. That is to say we assume the forward propagating fields are supplied, and allow the user to switch the direction if they want to. Use case: I record some fields in an output monitor, and I want to send these back in another simulation. I would just use the recorded FieldData and set direction = "-" in the CustomFieldSource.

Honestly though I'm also fine with us building the minimum viable product for now - just inject exactly what is given (tangential E and/or H), while other things like doing a numerical computation on H from E, or reversing the direction, can be done on the FieldData directly. I think I even prefer this now that I think about it...

[tylerflex]

Ok, next week I can set this up:

• remove CustomCurrentSource (?)
• usefield_data : FieldData as a Field()?

or how do you want to proceed? Note that we can still use .json files to test things before .hdf5 is supported, it just won't handle very large custom sources.

[tylerflex]

My last commit changed things according to our discussions:

• Removed CustomCurrentSource
• Added field_data : FieldData as the data container in CustomFieldSource.
• Added validator assuring a single frequency value in field_data within frequency_range of the SourceTime.
• Added validator assuring at least one tangential field component is defined (although maybe we dont need this as if not given, just inject nothing?)

A few things to do still:

• do we make monitor optional in MonitorData? If so, we need to change a few workings which assume that monitor is present. For example, discretizing for symmetry assumes self.monitor. Monitor is definitely a bit clunky in this CustomFieldSource object but useful in other places, so not really sure how to best handle this. I'm leaning towards making it optional and just making sure everything works without it as well as with it.
• Need to change apply_symmetry() to not take Simulation on both front end and backend.
• How do we want to handle the CustomMedium? what is the corresponding data object? For example, should PermittivityData contain sigma_xx and be the data container for the CustomMedium?
• backend implementation (of course)

[momchil-flex]

Overall, sounds good. Note that I already reverted apply_symmetry in my normalization PR #455. I still need to make this work on the backend, I was waiting for us to wrap up frontend discussions first. The only place where self.monitor is now used in a FieldData is in self.field_components. I wonder how we can replace this to automatically check what data attributes are present and are not None? Right now we don't have a list of the "possible" field names in a given data structure, although things like "Ex, Ey, Ez, Hx, Hy, Hz" are hard-coded in some places (symmetry values, grid locations...). Should we introduce a variable storing a list of the possible field component keys?

Now that we'll be using a FieldData for a source though, it definitely feels a bit clunky to have a monitor in there. I forget, is there any purpose of that - assuming we get rid of all self.monitor within FieldData itself? It seems to also be used just once (replacably) in load.py. Logically I wouldn't mind if FieldData is more purely a FieldData container, but I'm not really against keeping monitor but making it optional.

Regarding medium, PermittivityData can be complex, so I think we have everything we need? We just need to do the same thing for now: ensure that only one frequency is given. Then we can use that frequency to covert Im(eps) to sigma.

[tylerflex]

The only place where self.monitor is now used in a FieldData is in self.field_components. I wonder how we can replace this to automatically check what data attributes are present and are not None? Right now we don't have a list of the "possible" field names in a given data structure, although things like "Ex, Ey, Ez, Hx, Hy, Hz" are hard-coded in some places (symmetry values, grid locations...). Should we introduce a variable storing a list of the possible field component keys?

I think we can just hard code it in field_components() like we do in permittivity data.

class ElectromagneticFieldData:

    @property
    def field_components(self):
        all_components = dict(Ex=self.Ex, Ey=self.Ey, Ez=self.Ez, Hx=self.Hx, Hy=self.Hy, Ez=self.Hz)
        return {k:v for k,v in all_compoents.items() if v is not None}

Now that we'll be using a FieldData for a source though, it definitely feels a bit clunky to have a monitor in there. I forget, is there any purpose of that - assuming we get rid of all self.monitor within FieldData itself? It seems to also be used just once (replacably) in load.py. Logically I wouldn't mind if FieldData is more purely a FieldData container, but I'm not really against keeping monitor but making it optional.

Yea I dont know, the monitor type info is used to make the DATA_TYPE_MAP but we could also hard code that.. Also note that in SimulationData, if we know the monitor name (for accessing the data) we can also grab the monitor via sim_data.simulation.get_monitor_by_name(). I guess let's see how much of the need for monitor we can get rid of and go from there. maybe I'll try to rebase against develop and get rid of the simulation.discretize(self.monitor) dependence in the symmetry bits too.

Regarding medium, PermittivityData can be complex, so I think we have everything we need? We just need to do the same thing for now: ensure that only one frequency is given. Then we can use that frequency to covert Im(eps) to sigma.

Yea that should be fine, it's just slightly weird that regular medium requires real permittivity and a real conductivity, but I guess the permittivity data comes back as complex anyway so we might as well just go with that. should we try to do it in the same PR then?

[momchil-flex]

Yea that should be fine, it's just slightly weird that regular medium requires real permittivity and a real conductivity, but I guess the permittivity data comes back as complex anyway so we might as well just go with that. should we try to do it in the same PR then?

I think it would make sense to have some classmethods to define (n, k, f) or (eps, sigma, (f? because it's needed in the data)). If we only accept PermittivityData, the user needs to define eps_xx, eps_yy, eps_zz separately, which is tedious for most cases.

Probably separate PR makes more sense.

[tylerflex]

Ok, so I just made a few changes:

• monitors are optional in the MonitorData, reshuffled things around to avoid needing monitors as much as possible. will wait until we merge the normalize_index stuff in first to completely separate it.
• added custom medium and a few class methods, same PR just figured it would be simple enough to do here.

[tylerflex]

• Rebased against develop, fixed merge conflicts
• Fixed tests, note that I underscored _test_data_performance, it doesn't test much, breaks due to memory requirements of the GitHub actions and is more for internal benchmarking anyway.
• Made MonitorData.monitor optional. In SimulationData.apply_symmetry added condition to just return a copy if monitor is None (no way to discretize). Maybe we should handle this differently, I tried getting the monitor via monitor = self.simulation.get_monitor_by_name(monitor_data_name), but in some tests, the monitor was not in the Simulation.
• Added tests and imports for the CustomMedium

[momchil-flex]

• Made MonitorData.monitor optional. In SimulationData.apply_symmetry added condition to just return a copy if monitor is None (no way to discretize). Maybe we should handle this differently, I tried getting the monitor via monitor = self.simulation.get_monitor_by_name(monitor_data_name), but in some tests, the monitor was not in the Simulation.

So that's just because the tests set things up "wrong" no? The SimulationData is wrongly constructed for its associated Simulation?

[tylerflex]

So that's just because the tests set things up "wrong" no? The SimulationData is wrongly constructed for its associated Simulation?

Yea it was set up incorrectly (by accident) in test_sim_data:test_empty_io where the monitor name did not match the key in the monitor_data dict. To rectify this, I pushed a commit that does the following:

• added SimulationData validators that check whether each monitor_data element (a) has a monitor defined (b) is also in the simulation (c) has the same key as the monitor name.
• raises an error in apply_symmetry if the monitor_data.monitor is None (should already be checked by validator, but just in case). apply_symmetry works the same as before (accepts monitor_data, discretizes its .monitor) but only after the check.

This should be a bit more airtight now and we should be able to safely not require monitors in the MonitorData for "casual use" (if not used in SimulationData)