Screening currents near device boundary

[ark286]

Hi Logan,
I am working with Katja Nowack, and we are trying to use your superscreen package to simulate the diamagnetic response (e.g. susceptibility imaging) of some superconducting devices that we have measured with our scanning SQUID setup. I am having some issues with what I believe is some artifact in the simulation that is somehow related to the meshing. I'll try to outline the issue with a simplified example below.

Starting with the case where everything looks fine, a simple rectangle:

• I generate a landscape for lambda (x,y), which for the example here is assumed to be uniform

• I generate create polygons to match the geometry and mesh it using
  device.make_mesh(max_edge_length = 2,preserve_boundary=False,buffer=0)

• I model the FC current as a simple circle loop

• The FC current generates a magnetic field at the sample

• The field from the FC generates screening current in the sample

zooming in

• The screening currents in the sample produce a magnetic field at the SQUID plane that is threaded through the pickup loop, and we measure the flux

• I then raster the position of the FC current/pickup loop to simulate a scanning SQUID susceptibility image. Here, I am showing a horizontal line cut through the device.

Barring signal magnitudes, which I am not too concerned with at the moment, this all looks fine. Now if I alter the geometry slightly to begin to include leads for the device, I get some strange behavior.

Using a similar simple example:

-The lambda landscape (uniform in this case) and mesh now look like

• The FC current and position stays the same, but now the screening currents look like this:

zooming in

Notice the scale of the colorbar here as compared to before. There is now some additional current magnitude measured near the device edge.

• This results in additional field calculated at the SQUID plane from the screening currents and therefore additional flux through the pickup loop, again notice the scale of the colorbar here.

-If I simulate a similar horizontal linecut of the SQUID moving across the device, I now see an 'overshoot' of the susceptibility near the edges of the device.

This doesn't seem right. The addition of the leads shouldn't change the diamagnetic response of the device far from them.

The overshoot evolves with the horizontal length of the leads. As the leads get shorter, the overshoot tends to decrease and the left and right overshoot become more asymmetric. However, there is some 'turn on' behavior of this overshoot with respect to the leads i.e. if the lead is absent there is no overshoot, but if I make the leads one pixel long, an overshoot appears. My best guess is this has something to do with the meshing and perhaps something to do with the mesh points lying on the boundary i.e. maybe how boundary conditions (current conservation) is implemented?

I have tried a bunch of different things here, but the one thing that seems to remedy the overshoot is adding a buffer to the mesh:
device.make_mesh(max_edge_length = 2,preserve_boundary=False,buffer=10)

The linecut then looks much better, without the overshoot:

However, looking at the current density, there is now current flowing outside of the device boundary (ends at x = 5um) which is not physical.

Normally, I would just think the overshoot is some artifact of the simulation, which I still think it is, but we have susceptibility data measured on a device that shows a similar overshoot in the diamagnetic response. I think the overshoot in the data is stemming from some physics going on, but seeing the feature in the superscreen simulations made me think maybe there is a possibility of the overshoot being some SQUID imaging artifact. I realize this is a very long winded question, but if you have any insight here, it would be greatly appreciated.

Here's a link to a notebook with the example from above if that is helpful:
https://drive.google.com/file/d/1gu4jY5RsFgwSW8UP9ZQGX7nN6cKKA9ZX/view?usp=sharing

Cheers,
Austin

[loganbvh]

Hi Austin,

I played around with your model a bit. As far as I can tell, this is not a meshing artifact. Significantly increasing the mesh density near the sample edges (e.g. final_shape = final_shape.resample(801)) or over the whole sample (max_edge_length = 1) doesn't make the undershoot go away.

I think the undershoot is an effect of the device geometry combined with the very short effective penetration depth \Lambda ~72 nm. Weird things can happen near sample edges when the penetration depth is very short because the screening currents can be large and highly localized. The undershoot does get smaller if you make Lambda larger.

On the other hand, I agree that it's surprising that a small change in the sample geometry far from the field coil location can have a relative big impact on the susceptibility. It seems that the spatial distribution of the screening current is modified slightly by the presence of the leads. If you see a similar undershoot in the actual measurement, then maybe it is a real effect.

Some general suggestions:

1. Ideally the max edge length in the mesh should be smaller than the z-standoff at which you want to evaluate the field produced by the screening current (z_0 in your notebook).
2. Near sample edges, you want an even denser mesh, ideally with max edge length comparable to Lambda (although this is not always practical). For this, you can use Polygon.resample: polygon = polygon.resample(some_large_number_of_points)
3. While I provide the option to define a spatially varying Lambda(x, y), please be aware that (to the best of my knowledge) the London model is only formally correct and self-consistent when Lambda is constant as a function of position.
4. If you really care about the current distribution near the film edge (as you do in this case), then making the mesh with a buffer region is not a good idea. The reason is that the stream function g is defined at the mesh vertices. The current density is related to the gradient of g. The most natural way to define the gradient results in the gradient being defined at the center of each triangle in the mesh (as opposed to the actual mesh vertices). To evaluate the gradient at the mesh vertices, you take an average over all the adjacent triangles (see here). So when evaluating the current density, some current can "bleed" beyond the film boundary if there are any mesh vertices outside the film. If you had a much denser mesh near the edges, then this wouldn't be as much of a problem.
5. There is also just a limitation in terms of plotting the current density. When using matplotlib.pyplot.triplot, the color of each triangle is determined by the average value of the array over the three triangle vertices. This also leads to the appearance of current density outside the film boundary if there are mesh vertices outside the film. But it is just a plotting artifact.
6. For the sake of simplicity, you could try using a simple 1D current loop as the current source (example here). Also an easy way to make a pickup loop is to use a sc.Polygon: