s-wave collisions

[ElliotB256]

To simulate evaporation in optical/magnetic traps, AtomECS requires a way to calculate the effects of short-ranged (eg, s-wave) collisions.
Given that the collisions are short-ranged, all methods used should a spatial binning method for speed.
Beyond that, there are subtle variations in the implementation details, and we could support a number of approaches:

1. Monte-Carlo method: determine a number of collisions to occur within a box, eg based on density and mean velocity. Perform collisions between velocities of random atomic pairs.

   • How to do mixtures with this, eg fermions+bosons? Sr87+88?
   • Will be faster, because it's one loop per box with comparatively few collisions. Scales as O(N), where N is average atom number per box.

2. 'Pair-wise' Monte-Carlo method; Iterate over each pair within a box. Calculate the chance for a collision to occur, then draw a random number to determine if it did.

   • Easier to do mixtures, because you can specify collision lengths between pairs_
   • Will be slower, and scales O(N^2) where N is average atom number per box.

All models should also have the ability to specify some kind of 'macro atom', so that the behaviour of a large number of atoms can be represented by a numerically tractable ensemble. There's a few ways we could do this:

1. Have a simulation resource which specifies a global scaling parameter.
2. Specify per-atom MacroAtom components.

I like the global scaling parameter because it is straight forward. Also, I can see how we might be able to simulate multiple orders of magnitude of atom loss by rescaling the parameter as atoms are lost. For instance, this could allow one to start with 1e4 atoms representing 1e10 atoms and end with 1e4 atoms representing 1e4 atoms.

[ElliotB256]

Worth noting that collision length isn't really a per-atom quantity - it's defined between a pair of atoms - so we can't really have a component for it.

How to store collision lengths for calculations? I would suggest we have a table which holds the collision lengths for all species, and index this table using per-atom collisions::Species { id: byte }. The table should form part of the Resource option which is used to enable/disable s-wave collisions.

[ElliotB256]

@d-garrick Note that a more general implementation can't just swap velocities, it needs to swap momenta in the COM frame.

Also, it can't just swap momenta in the COM frame either - otherwise that would imply all x-components of velocities are decoupled from all the y-components, which is unphysical. Look up references and prev work to see what model is best to use.

[MauriceZeuner]

Hi @d-garrick ,

I'm Maurice and I'm currently doing the dipole-force branch. Since we are ultimately interested in simulating evaporation (for example in a dipole trap) - I wanted to ask if you are interested in a brief meeting about the collisions. I have implemented a really lame version of them in python beforehand (just to get a feeling for it) and I'd be interested to learn how you implemented them and how to use them.
Of course, I can also give you the outline of how the dipole force will work - although it is very similar to the cooling light and forces functionality already in the master branch.

Anyway, I think it would be great to meet so we are aware of what everyone is doing, what do you think?

Greetings,
Maurice

[d-garrick]

Hi Maurice, We can have a meeting about this - the collisions code is still unfinished and requires more testing, but I'm happy to have a chat about it with you. I can probably chat on teams tomorrow (Thursday) or Friday. When are you free? David

…

------ Original Message ------ From: "Maurice Zeuner" ***@***.***> To: "TeamAtomECS/AtomECS" ***@***.***> Cc: "d-garrick" ***@***.***>; "Mention" ***@***.***> Sent: Tuesday, 23 Mar, 2021 At 20:52 Subject: Re: [TeamAtomECS/AtomECS] s-wave collisions (#14) Hi @d-garrick <https://github.com/d-garrick> , I'm Maurice and I'm currently doing the dipole-force branch. Since we are ultimately interested in simulating evaporation (for example in a dipole trap) - I wanted to ask if you are interested in a brief meeting about the collisions. I have implemented a really lame version of them in python beforehand (just to get a feeling for it) and I'd be interested to learn how you implemented them and how to use them. Of course, I can also give you the outline of how the dipole force will work - although it is very similar to the cooling light and forces functionality already in the master branch. Anyway, I think it would be great to meet so we are aware of what everyone is doing, what do you think? Greetings, Maurice — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#14 (comment)> , or unsubscribe <https://github.com/notifications/unsubscribe-auth/AST4AXYLLOK6NJQGS55WAVLTFD5QHANCNFSM4XY6IOCA> .

[MauriceZeuner]

Hi David,

great, thanks! I am very flexible tomorrow - Friday before noon as well. Would you like to suggest a time that works best for you?

Best wishes,
Maurice

[d-garrick]

Shall we say 11am tomorrow on Teams? David

…

------ Original Message ------ From: "Maurice Zeuner" ***@***.***> To: "TeamAtomECS/AtomECS" ***@***.***> Cc: "d-garrick" ***@***.***>; "Mention" ***@***.***> Sent: Wednesday, 24 Mar, 2021 At 19:30 Subject: Re: [TeamAtomECS/AtomECS] s-wave collisions (#14) Hi David, great, thanks! I am very flexible tomorrow - Friday before noon as well. Would you like to suggest a time that works best for you? Best wishes, Maurice — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#14 (comment)> , or unsubscribe <https://github.com/notifications/unsubscribe-auth/AST4AX5F7KXNLWTHEE7L6QLTFI4TTANCNFSM4XY6IOCA> .

[MauriceZeuner]

Sounds good - I think we don't have a common teams channel yet. You should be able to find me by my name "Maurice Zeuner" or ----.

See you,
Maurice

[ElliotB256]

This graph shows agreement of the MC approach with theoretical rates as a function of collisions/timestep/box and atoms/box. Note that even 30% difference is expected due to inaccuracies in the numerical calculation of the 'expected' rate (I'm using only approximate factors from kinetic theory, there's a few factors missing).

[ElliotB256]

(The above simulation is for a box trap simulation domain with 10 boxes along each axis for 1000 boxes total.)

[ElliotB256]

Slightly more correct version of the same:

[ElliotB256]

@MauriceZeuner

[MauriceZeuner]

Is there a minimum example that shows how to use the collisions? I'd be very interested to throw some particles into a dipole trap and evaporate.

[ElliotB256]

@DavidGarrick is still working on it

[d-garrick]

I can add a simple example to the collisions branch tomorrow. It produces accurate collision rates for thermal clouds but haven't tested evaporation yet.

[ElliotB256]

@DavidGarrick I know we spoke a few times today, could you post a brief update here? :)

[d-garrick]

I've been testing collision rates and evaporative cooling and had some discrepancies but have found the bug - the current version overstimates the collision rate by a factor of the macroparticle value. After fixing this, the s-wave collisions (and magnetic trapping) seem to work well and agree qualitatively with experimental data:

This is a phase space density curve of evaporative cooling in the Time-Orbiting Potential trap, which I've also implemented on the magnetic_force branch. This experimental data is quite old and not perfectly suited to comparison to simulation, so I'm hoping to take some fresh data explicitly for the purposes of comparing to simulation, but the agreement is still reasonably good.

The model for collisions we use assumes that we have identical particles, that the velocity distribution in the collision cell is at least nearly thermal, and that we can model a cell as a uniform density. For most scenarios and a suitably small cell width these last two are probably good assumptions. We may want to modify the existing code to handle mixtures of species.

I also noticed some other interesting bugs. One is that the collisions must be carried out after the velocity verlet integrator has finished, or else the integrator will no longer conserve energy (in quite a major way).

[d-garrick]

It's also worth mentioning that at the moment at least, the user has to spend some time figuring out good parameters to use to ensure accuracy of the simulation. For example, the timestep must be significantly less than the average (or possibly the smallest, in the case of inhomogeneous density) collision time, at least 10 times less. I have also read in the literature that the timestep needs to be much smaller than the mean transit time; that is, the average time taken for a particle to cross a collision box. However, I have not personally noticed any unphysical effects due to this, and I am actually not sure why it would matter on a theoretical basis.
The collision box width needs to be large enough to have at least >10 but preferably >20 particles, but too large and the simulation will be slower and correlations between position and velocity will be increasingly coarse-grained. Angular momentum is not conserved in DSMC collision algorithms, so any simulation in which this is of interest is likely to be inaccurate.
Small particle number per box is of course fast, but leads to large statistical fluctuations. A time-averaged (or spatially averaged, over many identical collision boxes) collision rate will be accurate, but the fluctuations need to be small enough to average to close to the correct value faster than any relevant dynamics in the trap.

[ElliotB256]

the timestep must be significantly less than the average ... collision time, at least 10 times less.

Can you produce a 2D plot like the one above to show when it breaks down?

Graph looks great!

[ElliotB256]

@DavidGarrick discussed in slack the advantages of making a PR for these early collision features, with more configurations exposed/optimisations in a second PR.

[ElliotB256]

@DavidGarrick ready to pull?

[MauriceZeuner]

Hi David,

I think you mentioned some unit-tests still failing - just to check:

collision_rate() --> fails
test_collisions() --> deaktivated/fails if ran

Also, there is a panick that I don't understand where it actually occurs:

<unnamed>' panicked at 'okTried to fetch resource of type 
`CollisionsTracker`[^1] from the `World`, but the resource does not exist.

Is that the currently "correct" behaviour or did I mess something up while merging?
Thanks :)

[d-garrick]

The unit tests are very unfinished so I can't comment much on which particular ones are failing, they won't be done until I try and merge collisions (which I'll try and work on soon in between lab work).
Let's email about any errors or questions and if we find an actual bug we can mention it here.

[ElliotB256]

I saw @d-garrick working on this again today

[ElliotB256]

Just out of interest (and when the experiment is operational again) - did you ever take a more detailed set of data for the evap ramp and try fitting to see what the most realistic cross section for Rb is? eg, ~100 different values for atom number, find the cross section which gives best agreement with experimental data. Then compare that to the expected cross section from theory?

[d-garrick]

I haven't taken that data - might take it once the experiment is working again. Could try, rather than evap, damping rate of a breathing oscillation?

[ElliotB256]

Any chance to revisit that data now the experiment is operational?

[ElliotB256]

Tagging @BrianBostwick (hope it lets me)