Dipole Trap asymetrical heating

[MauriceZeuner]

Hi all,

sorry for being quiet - I have been desperately trying to find the cause for a heating effect which seems to be still in the code but is most likely a dipole trap specific feature. I apologize in advance, this is super nasty and I'd be glad for new ideas how this could happen.

First observation: If you simply put atoms in a dipole trap, but also have some CoolingLight Entities in the game, however, very far detuned, so they shouldn't do anything - you'll notice that there appears to be a preferred direction for their escape - the x-direction. It can best be observed in the (branch: new-dipole-force) cross_beam_dipole_trap.rs example - where there should be no difference between x and y since we have two identical, horizontal beams in exactly these two directions (and the MOT beams are on the diagonals). However, if you plot the absolute coordinate values of the ensemble (or just have a very sharp eye when plotting the data in 3D in Blender) you'll see that they are more likely to be further away from the centre in (plus or minus) x-direction.

In addition to that, there is still a good amount of numerical heating happening on top of that, this effect, however can be reduced by increasing the detuning of the MOT beams even further:

Before you ask: I checked in a pure symmetrical 3D MOT if there was a difference between x and y and I don't think it is present there. Of course, I also checked the way the dipole force is calculated for any symmetry breaking but couldn't find one. At that point, however, I am running a bit out of ideas where to look for it.

For comparison, here are the same diagrams as above just without any CoolingLight entities in the code:

How did I find this behaviour in the first place? Well, I was running MOT -> DipoleTrap transitions with three stages:

1. MOT only (first 0.2 s)
2. MOT ramp down (power and detuning) (another 0.2s)
3. Dipole only - CoolingLight entities are deleted - (rest 1.4 s)

(see file: example/red_mod_xodt_transition.rs)

which looks like this: and

I'd be extremely grateful for any, even far-fetched, ideas what could cause this behaviour.

[ElliotB256]

however, very far detuned, so they shouldn't do anything

what detuning are you using? perhaps you could quick plot the offset versus detuning?

[ElliotB256]

Also, to verify that the cooling light is the cause - what happens if you change the direction of the cooling light?

[MauriceZeuner]

what would you change it to?

[MauriceZeuner]

Hang on. Sorry, this does exist on the master branch and for the red MOT also:

[MauriceZeuner]

See branch asymmetric-heating and then example file: 3d_red_mot_from_center.rs

this is the master branch, I just added that example file to it 36da241

[MauriceZeuner]

Ha! But not for the blue MOT!

So it's a red MOT specific problem.

[ElliotB256]

Can you provide a minimum example that produces the behaviour?

[MauriceZeuner]

yes, already have: 36da241

[MauriceZeuner]

plug in the transitions for the Blue MOT, change power and detuning to something that makes sense and you don't get the effect. It's only for the red MOT, it seems. If you don't trust the CentralCreator, I can also explicitly create random atoms in the center, but I think that's not the reason. (Since there's no problem with the blue MOT)

[MauriceZeuner]

Found, that the x-direction receives about 26% more force than y. By that, I mean the summed, absolute values in each direction. In a symmetrical setup, the forces should be equally strong. I investigated it using: a1c0f12

[MauriceZeuner]

I think I'm getting closer. RateCoefficients seem to be imbalanced too, not by much, though. However, something really strange happens in z-direction:

WITHOUT gravity ┌ ┐ │ 4396097345004.737 │ +z │ 3880732997761.6094 │ -z │ 4601255035803.305 │ +x │ 4631279125587.867 │ -x │ 4417368384127.704 │ +y │ 4331600239418.176 │ -y └ ┘

WITH gravity:

┌ ┐ │ 5236480448946.082 │ │ 2962149027246.679 │ │ 4526613812299.344 │ │ 4580436843213.253 │ │ 4272408756990.435 │ │ 4307329865720.1665 │ └ ┘ see 4b53446

II will investigate further.

[ElliotB256]

It might help to measure the rate coefficients for atoms placed at well-defined positions?

[MauriceZeuner]

Sorry, have to retract my last observation. Still had a 100µm offset in z for the magnetic field. RateCoefficients are sound (just a few examples):

pos: Vector3::new(5.0e-4, 5.0e-4, 5.0e-4) ┌ ┐ │ 7835.07551271332 │ │ 24398.235402460607 │ │ 16297.298995739506 │ │ 24578.87894061315 │ │ 16297.298995739506 │ │ 24578.87894061315 │ └ ┘

pos: Vector3::new(5.0e-4, 0.0e-4, 0.0e-4) ┌ ┐ │ 486199.11196358874 │ │ 486199.11196358874 │ │ 9044.374930279317 │ │ 1911505.4929097523 │ │ 486199.11196358874 │ │ 486199.11196358874 │ └ ┘

pos: Vector3::new(0.0e-4, 0.0e-4, 0.0e-4) ┌ ┐ │ 31767.037170127885 │ │ 31767.037170127885 │ │ 31767.037170127885 │ │ 31767.037170127885 │ │ 31767.037170127885 │ │ 31767.037170127885 │ └ ┘

But the heating and asymmetry remain.

[ElliotB256]

can you comment on whether the numbers above are expected? particular element 0 of the first one, and elements 3 and 4 of the second one

[MauriceZeuner]

yes, everything is as one would expect. the position is given above the vectors and the elements denote the entries of the RateCoefficients component. The first two are plus and minus z, 3 and 4 are x, 5 and 6 are y

[MauriceZeuner]

Just placed a single atom in the trap centre (static!). Summed force over 100000 time steps is:

┌ ┐ │ 0.0000000000000000014452095452848004 │ │ 0.0000000000000000011562930341955761 │ │ 0.0000000000000000010981273888407856 │ └ ┘

meaning, x is still getting more force even in that example. So we know it is somewhere between the RateCoefficients (which work fine) and the Force - which does not work fine. Also, it is (at least not alone) the integration system.

[MauriceZeuner]

It's not the ExpectedPhotonsScatteredVector either, that one is fine!

That only leaves the drawing from the Poisson distribution and the two ForceSystems. I'll look at them tomorrow.

[ElliotB256]

Looking at these vectors with all that wasted precision does make me feel sad :D I wish the uom crate was able to work with Vectors. Maybe I can do something with macros to get it working - it's a really nice crate if we can use it, and it also does compile-time dimensional analysis.

[MauriceZeuner]

What wasted precision do you mean? The leading zeroes? It's just a println! command, not a value from a debugger (VS Code's debugger does not work for some reason...)

[MauriceZeuner]

Okay, it appears to be a problem that exists when both (!) the EmissionForceOption and the ScatteringFluctuationsOption are active. Summed force vector for static, central atom after 10,000,000 steps:

no emission force, only scattering fluctuations ┌ ┐ │ 0.0000000000000000732180549846602 │ x │ 0.00000000000000007321997707152613 │ y │ 0.00000000000000007307005429598201 │ z └ ┘

only emission force, no fluctuations ┌ ┐ │ 0 │x │ 0 │y │ 0 │z └ ┘

both! ┌ ┐ │ 0.00000000000000016649020096496058 │x │ 0.00000000000000013291046326512055 │y │ 0.00000000000000013254777955999218 │z └ ┘

We can see that only in the last case (both activated) the x-direction gains significantly more force.

[MauriceZeuner]

Found the (or at least one) guilty! Its maths::random_direction(). For some reason it is more likely to spit out a vector with a larger absolute value in x.

[MauriceZeuner]

And the reason why it is important in red but not blue is that in blue, the function is called rarely since there are usually more than 5 photons scattered per step.

[ElliotB256]

ahh, interesting! is there a bug report for random_direction you can link?

[ElliotB256]

(and are you saying the discrepancies go away if you turn off fluctuations?)

[MauriceZeuner]

Yes, if you simply set laser::force::EmissionForceConfiguration::explicit_threshold to 0 (for the EmissionForceOption), the asymmetry is gone, hurray! Not sure if all problems in the dipole trap are gone but in the master, that should solve the problem. To solve it a bit more sustainably, I recommend using a different method to obtain random vectors. Looking at the function it is actually quite obvious why it does not produce homogenous results.

What do you mean by bug report? Is that another github feature I don't know about?

[ElliotB256]

Hmm, I didn't realise we were defining our own function for it! I assumed we were using one from another crate. Can you replace ours with https://docs.rs/rand_distr/0.4.0/rand_distr/struct.UnitSphere.html ?

[MauriceZeuner]

yes, will do that. Should I make that change directly on master?

[MauriceZeuner]

Hm, the heating for the dipole trap remains, but at least it is now symmetrical.

So I guess, problem 1 out of 2 solved.

[ElliotB256]

If it only occurs when you do actually have random kicks drawn (ie, photons scattered), then are you sure the heating is unphysical? There is still some scattering rate for a far-detuned trap - what is the scattering rate you expect given your parameters and what is the associated heating rate?

[MauriceZeuner]

yes, because I did not give my dipole Lasers a "CoolingLight" component so far - so they should not be able to scatter at all, I think.

[ElliotB256]

Can you verify if it is integration error by reducing timestep by a factor of ten and running the same simulation time?

[MauriceZeuner]

I tried what you suggested and it makes no notable difference. It's not an integration error, I'm pretty convinced. Just for reference:

A: startup phase, red MOT only B: Switch on Dipole, start ramping down MOT C: Dipole only, switch off red MOT completely

diagram shows distance from center for all particles

[MauriceZeuner]

My best guess is that something goes wrong when I switch off the MOT

[ElliotB256]

Sounds like you'll have to go through the systems one by one until you find the cause

[MauriceZeuner]

Found the problem: In short, the problem was that various components attached to the atoms were not cleared when the CoolingLight components are deleted. The ApplyEmissionForceSystem now runs despite no CoolingLight components (entities with them) are present anymore because it takes that information from the ActualPhotonsScatteredVector which still contains non-cero values from the time when the Laser beams were still active.

Ultimately I propose, we should create a clean, general System or method that clears component data of entities that are deleted, especially beams. And, for runtime reasons, also leads to irrelevant systems not running without any use.

So far I have only made the hotfix to my dipole example code which basically turns off the ApplyEmissionForceSystem after the MOT-beams are deleted.

[ElliotB256]

Good find! I think a few things are confused here though.

First, if you delete an entity the associated components are removed once you call world.maintain. For example, from the delete doc:

Deletes an entity atomically. The associated components will be deleted as soon as you call World::maintain.

It sounds like the bug is really in the formulation of the system logic. Looking at the systems, the following system clear out the values at the start of the frame:

pub struct InitialiseExpectedPhotonsScatteredVectorSystem;
impl<'a> System<'a> for InitialiseExpectedPhotonsScatteredVectorSystem {
    type SystemData = (WriteStorage<'a, ExpectedPhotonsScatteredVector>,);
    fn run(&mut self, (mut expected_photons,): Self::SystemData) {
        use rayon::prelude::*;
        use specs::ParJoin;

        (&mut expected_photons).par_join().for_each(|mut expected| {
            expected.contents =
                [ExpectedPhotonsScattered::default(); crate::laser::COOLING_BEAM_LIMIT];
        });
    }
}

If no lasers exist, this values should remain at the default (zero) value. However, it seems this system wasn't added to the dispatcher! Could you check into this for master?

[ElliotB256]

(It's fully possible I removed the system by accident during optimisation aswell! :D )

[ElliotB256]

Closed because bug was found and new issue opened