Replicate single particle electron motion from Gummersall

[MatthewGrim]

This issue aims to replicate some of the single particle motion results from a paper by Gummersall using my previously implemented PIC methods.

[MatthewGrim]

To do this, I first need to carry out some grid convergence studies. The following replicates the Gummersal conclusion that $\Delta t$ should be $10^-9 \times radius$.

Study Results

Me playing around

The first of these is my own study playing around with dt so that it was some factor of the time taken to traverse the device.

This shows that there is little to no convergence at these radii.

Replicating Gummersall

The second study considers dt factors relative to the Gummersall hypothesis. Here there is much better convergence.

[MatthewGrim]

Field Convergence Study

Completing a simple field convergence study, it seems that 200 is a relatively good resolution to choose. There are regions with significant deviation, but I assume these are close to the loops, where particles are not going to reach anyway.

Update plot

[MatthewGrim]

Plan for the weekend

Now the issue I had with the interpolator is fixed, I can try to generate meshes. I have resolved that errors in the B field approximately fade after 200 points. So I need to generate these meshes for future studies. For my simulations, I will consider only 10cm, and 1m radius devices. Larger devices are not reasonable for prototypes of polywells so for now to keep the study as lightweight as possible, I will try to generate these meshes for 100, 1000, and 10000 Amps. 6 meshes for future studies. The meshes at 10cm will be used to replicate figure 2 from Gummersall.

So the plan is to build the 6 meshes over the weekend. Carry out a simple grid covergence study with the resulting mesh. Then replicate the result from Gummersall over the next week.

I wanted a domain resolution of 200 points, but this is too difficult to generate over a weekend in python. Instead I am going for 130 points, arbitrarily.

In the future

I should consider speeding up simulations by making the system periodic, and writing a Fortran version.

To be less data intensive, I should convert the text file outputs of the vectors into a vtk array.

[MatthewGrim]

Results of grid convergence with mesh 130 points in domain resolution

The results are similar to previous studies. Convergence is achieved for the 10cm radius mesh, sometimes at earlier an order of magnitude earlier than 10^-9. This applies for the 100A, 1kA, and 10kA cases.

High velocity simulations

Going to the high end of the electron velocity side, it is hard to get convergence between the simulations for high currents. The particle is always escaping, but the positions generated vary. Lower resolution simulations seem jagged, so they are definitely not suitable.

In terms of position, at low current density, the position of the lost electron is similar. At high current density, the particles leave out of different cusps. This indicates I may need to add adaptive time stepping to keep the simulation resolved when moving through the cusps, where there are very tight oscillations at high speed.

[MatthewGrim]

First pass at replicating result

This first pass was done without varying the initial particle position from the origin. The time scales are slightly longer than in the plot by Gummersall.

[MatthewGrim]

Digitising the plots from Gummersall, I see quite a significant discrepancy:

The reason for this is likely to be because of the position of initialisation for the electrons. Being positioned at the origin, it's likely that on average it will take slightly longer for the electrons to escape the device, because they have to drift through the whole radius to get out of the device. I will repeat the study with the initialisation similar to Gummersall and see how that changes things.

With a random spacing in the distance 3a/16, I get:

[MatthewGrim]

Replication at Gummersall grid spacing

Gummersall says in his thesis that he uses a much lower resolution for his mesh - 20 points across the radius of the device. This could be the reason the perceived fields are stronger in his simulations? I'm going to check if the results match when I use this resolution.

[MatthewGrim]

Results are still off for lower resolution sims that match Gummersall's mesh resolution. What are the other differences between our approaches?

[MatthewGrim]

Comparison with results using 200 resolution current loops without interpolation

This is treated as the ground truth. The one that matches this will be the one which adds up to me.

Result

Plotting all of the test runs I have done so far at different resolutions:

I am converging to a different value...

OK - the next step is to try the particle pusher discussed in the paper. This pusher doesn't seem to be handling the gradB drift I would have expected in a Polywell device, which is why I would have kept to the standard Boris solver, like Gummersall seems to do in his thesis.

Before doing that i will regenerate the data to make sure I'm not being stupid...

[MatthewGrim]

Problem with uniform velocity and position distributions

I think there might be a problem with the way I am generating the velocities, they might be weighted in particular directions.

[MatthewGrim]

The uniformity of the velocity distribution is now fixed.

[MatthewGrim]

Even with the velocity being uniform there is still a significant discrepancy between our results. I am now going to try the simple particle pushed written in the paper.

Update There seem to be some discrepancies between the definitions in the paper, and in Gummersall's thesis. It's unclear whether he is using a Boris solver, or the analytic evolver, which I don't back. It also is unclear if there is any time varying B field in the plot I am trying to replicate....

[MatthewGrim]

Definition of escape was wrong

I just realised my definition of escape has been incorrect. It should be when the particle meets the loop planes. For simulations so far it has been when the particle leaves 1.1 times this value. This hopefully will correct the error in simulation results, but I have to wait until Monday to find out.

[MatthewGrim]

Run with escape definition changed

Even with the change to the domain size, the simulation results differ significantly. It's strange that the plots align so well with values one order of magnitude greater than them.

Update to plot including 1e5 curve

[MatthewGrim]

Replicating radial probability density results

In the last few weeks, I have been implementing storing histograms of the radial probability density distributions. These results look similar to what is below:

At low current, the electrons are escaping, but in the range between 10 and 100kA we get a better confinement.

In terms of the scaling of average radius, this is very similar between devices of different sizes, yielding average radius plots similar to what is below:

This differs considerably from the plots in Gummersall's paper again...

[MatthewGrim]

Conclusions

Qualitatively, the results I obtain are similar in many respects to Gummersall. At the same time, the results differ considerably in the final output - the radial probability distributions. It's not clear to me why this is, but I think I need to re-evaluate my simulation methods again, in order to make sure they are validated.

A silver bullet answer would be to implement the incorrect evolver seen in Gummersall's paper and see if this replicates his results. If it does, I am not concerned about the difference. If it doesn't, and matches mine ... I would start to be concerned about my mesh generation functions, and their accuracy ... despite having validated them to the best of my ability already...

I am closing this issue and opening two new ones - the first will re-evaluate the Boris solver, and mesh generation functions to ensure they are correct. The second will get more data out of the simulations I am running. I will try to get the velocity distributions, number of escaped particles, and their escape locations.

[MatthewGrim]

Additional Replication Plots

Replication of Gummersall Figure 5

Replication of Gummersall Figure 6

[MatthewGrim]

Comment on above results

Both figure 5 and 6 show variations from Gummersall's initial results that again make me doubt my methods, and mean I need to do some more digging. At the same time, the trends are similar (the gradient in figure 6 is extremely similar. Similar to figure 2, I wonder what would happen if I made the current 100kA instead of 10kA in figure 6?