atcollab / at

Accelerator Toolbox
Apache License 2.0
48 stars 31 forks source link

GWigSymplecticRadPass is wrong #190

Open lfarv opened 3 years ago

lfarv commented 3 years ago

The introduction of the the wiggler contribution in the computation of radiation integrals gives two independent methods to compute the damping times:

  1. atx derives the damping times directly from the 6-D one-turn transfer matrix (pure tracking)
  2. atsummary/ringpara use the radiation integrals. In the case of a simple planar, infinitely wide, single-harmonic wiggler located in a non-dispersive section, the computation is fully analytical.

For running the comparison, we used the lattice provided by @TeresiaOlsson here, and activated a single wiggler, the one labeled I12. @mashl showed here that to get a good description of the dispersion inside the wiggler, Nstep > 20 is necessary, so we tested Nstep = 5 (default), 20 and 40, without any significant difference. Here are the results:

Bare lattice

Damping times

x z s
atx 0.014234 0.019409 0.012004
atsummary 0.014213 0.019527 0.012008

Damping partition numbers

x z s
atx 1.370257 1.004900 1.624844
atsummary 1.373869 1.000000 1.626131

The bare lattice, without wigglers shows a very good agreement between the two methods:

With I12 wiggler

Damping times

x z s
atx 0.009844 0.019029 0.012336
atsummary 0.012210 0.015936 0.009403

Damping partition numbers

x z s
atx 1.727657 0.893717 1.378626
atsummary 1.305122 1.000000 1.694878

The results are here very different. In particular the vertical damping partition number from atx is quite different from 1.

An even simpler check can be done by comparing the synchronous phase from the 6th coordinate of the 6-D closed orbit (tracking) and from the analytical computation.

Bare lattice

synchronous phase [rd]
findorbit6 2.585603
atsummary 2.585573

Good agreement again.

With I12 wiggler

synchronous phase [rd]
findorbit6 2.585628
atsummary 2.438302

The synchronous phase from tracking does not change when introducing the wiggler, as if GWigSymplecticRadPass would not give any energy loss.

All these results indicate that there is something missing in GWigSymplecticRadPass.

lfarv commented 3 years ago

@willrogers: can you trace back where the GWigSymplecticRadPass integrator comes from? If we cannot get it repaired, I think we'll have to remove it from the distribution since it is obviously wrong… I can look inside to try and understand the problem but it may be long and unsuccessful !

willrogers commented 3 years ago

OK, I will see what I can find.

willrogers commented 3 years ago

It was initially added by @mashl here, although it was merged in a different pull request later.

mashl commented 3 years ago

@lfarv , @willrogers Hello, I find a bug on the gwigR.c code, maybe it could be helpful. The 324th line of the code should revised. current version : (line 319-324) if (fabs(ky/kw) > GWIG_EPS) { syky = sin(ky y)/ky; } else { syky = y sinc(ky y); aypx = aypx + (pWig->VCw[i])(kw/kz) kxshxsykysz; }

revised version: if (fabs(ky/kw) > GWIG_EPS) { syky = sin(ky y)/ky; } else { syky = y sinc(ky y); } aypx = aypx + (pWig->VCw[i])(kw/kz) kxshxsykysz;

and there is same problem at gwig.c

lfarv commented 3 years ago

@mashl : You are right, I'll make the tests again. Not sure it will solve the problem of missing energy loss, but it may answer another doubt I had concerning tune changes. Give me a few days to look at that…

lfarv commented 3 years ago

@mashl , @willrogers: So on bug less, but not the one I'm looking for. No change from the results I showed above: the synchronous phase derived from tracking is still identical with and without wiggler, meaning that there is no energy loss in the wiggler !!

mashl commented 3 years ago

dear @lfarv I see the change in vertical tune (phase advance) when the planar horizontal wiggler is located in ring, but the horizontal tune is same. do you expect the change in H-tune as well?

mashl commented 3 years ago

for instant: (all computations are performed by using ringpara) case 1 : planar horizontal wiggler ID : atwiggler('WIG', LWIG, LWIG/NP,B0, E0,'GWigSymplecticPass'); LWIG=2; B0=1.5; NP=40; Bare lattice tune: [0.1610 0.2247] Lattice with ID: [ 0.1610 0.2306]

case 2 : multi-harmonic elliptically polarized ID ID: previous ID with below modification ring1{118}.NHharm=2; ring1{118}.NVharm=2; ring1{118}.By=[0.7,0.7;0.8,0.2;0,0;1,2;1,2;0,0]; ring1{118}.Bx=[1,1;0.8,0.2;1,2;0,0;1,2;0,0];

Bare lattice tune: [0.1610 0.2247] Lattice with ID: [0.1818 0.2287]

mashl commented 3 years ago

Oh sorry I check betatron phase advance

lfarv commented 3 years ago

@mashl: in the case of a simple wiggler (planar, single harmonic), the tune shift can be analytically computed. For an "infinitely wide" wiggler (kx=0), there is only a vertical tune shift (what you noticed is perfectly correct). My problem is that the tune shift from tracking is almost twice the analytical value. But this is a minor problem compared to the missing energy loss. Now concerning the bug you fixed, it concerns only vertical (or elliptical) wigglers, which I did not look at yet.

willrogers commented 3 years ago

I just found a comment from Marco, one of our accelerator physicists, from July:

GWigSymplecticPass does not have the energy loss, no. well the story goes that we then found and tried GWigSymplecticRadPass, which seems to take into account the proper energy loss (good) However, when we check the beta functions, betax is devastated, not just a perturbation on the usual function. So we are not sure what to use and how,

  • without Rad we have a credible change in the vertical beta/ tune
  • with Rad we can see a good energy loss but the optics become unreliable, or so it seems

hence the idea of sharing experiences. I appreciate the difference in energy can put our machines on a different ground.

mashl commented 3 years ago

Dear @lfarv , Effects of ID- atx, ringpara, atsummary.docx you can find the comparison between atx , ringpara and atsummary for bare lattice and lattice with a horizontal simple wiggler in the attached doc. the results of atx and ringpara for all parameters seem reliable, but I think atsummary needs many modifications to be compatible with the presence of the ID. It should be noted I changed the ringpara and atx and they are different from the master branch codes, for instance, the diffusion matrix of a wiggler is added to the ohmienvelope code and findorbit6 is changed. SR_02_01.txt SR_01_01.txt

lfarv commented 3 years ago

Hello @mashl. Thanks for the full dataset, so that I could reproduce your tests ! You must be using an old version: since the merge of #188 (October 12), ringpara and atsummary share the same updated code and give exactly the same results, in agreement with your RINGPARA column, except for the natural emittance, where both ringpara at atsummary give 2.50681E-10 m. Concerning atx, I am interested in looking at what you did for the diffusion matrix. In which branch are you working?

mashl commented 3 years ago

hello @lfarv , I just update my branch (fork), and it's the same as what I use on my pc.