IdTablePass

[oscarxblanco]

Dear all,
dear @lfarv , @lcarver , @carmignani , @simoneliuzzo

This branch implements a new element, called InsertionDevice, in pyat. It uses the pass method IdTablePass for tracking, and requires : a name, the number of slices, an input file with a magnetic field map from Radia, and the beam energy in GeV. For example,

u20 = at.InsertionDevice('u20',10,"u20_g45mm_kicks_2022nov10.txt", 2.75)

The InsertionDevice element has been checked against AT MATLAB results for a particular magnetic field map, energy and ring, giving the same tune variation and frequency map.

 ... Tune variation ...
expected difference  : 0.0000 0.0044 -0.0000 (from matlab)
calculated difference:[ 7.00984399e-07  4.35046416e-03 -7.11769372e-10] (from python)

The InsertionDevice element, however, does not have yet a function to be read or written, from or to, .mat or .m files.
In addition, the function lattice.insert is incompatible with the InsertionDevice parameters because insert does not copy the Energy property of the element. As a workaround, I have created a dummy element in the lattice that is overwritten by the insertion device.

# create a new lattice where we will insert the InsertionDevice
newlattice = lattice.deepcopy()
### 2023feb05: insert function does not work, it removes the Energy property
### newlattice.insert(insertIDindex+1, u20)
### work around because at.lattice.insert does not copy the Energy property
dummyelem = newlattice[insertIDindex].deepcopy()
newlattice.insert(insertIDindex, dummyelem)
newlattice[insertIDindex] = u20.deepcopy()
# substract half length on each side
newlattice[insertIDindex - 1].Length = newlattice[insertIDindex - 1].Length - u20.Length/2;
newlattice[insertIDindex + 1].Length = newlattice[insertIDindex + 1].Length - u20.Length/2;

This pull request is related to issue #522, and the kickmap method mentioned in issue #1.

Best regards,
o

[lfarv]

Nice addition, thanks
Given the size of the code, I think you could move it to a dedicated module elements.py now really becomes very long…

[simoneliuzzo]

Dear @oscarxblanco,

yes, very nice addition indeed!

Yesterday a colleague from the ESRF ID group (Gael) asked me to check the impact of a wiggler on our lattice. Perfectly synchronized to your pull request!

I used your very useful demo lines of code in the header of this pull request to introduce the wiggler in my lattice model with very minor changes (my ID location is a marker).

I have used 3 radia kickmaps at different gap apertures for a wiggler. The kicks are order of 1e-9 rad when the wiggler is open and order of 1e-6 with the wiggler is closed. (see images)

What are the expected units for the kickmap? radians? I see that in the code the kicks are multiplied by 1/(Brho^2), may be this is how it is supposed to be.

What are the expected units for the positions? meters?

Could you share the kickmap file in your example?

I am trying to look at the results and I am still not really understanding them. I may have done some mistakes.

For example the optics are modified only around the wiggler (index 1806). The rest of the lattice is totally unaffected (0.0 % beta-beating and dispersion deviation).
Also, even if the kicks are changing by several orders of magnitude the impact on optics is (sharply) identical in the 3 cases.

The ID elements are filled, and the values do correspond to those in the input files (considering the 1/brho^2 factor):

print(w150_25mm)
InsertionDevice:
FamName : w150_25mm
Length : 1.5
PassMethod : IdTablePass
Filename_in : kickmap_w150_25mm.dat
Energy : 6.0
Nslice : 150
xkick : [[ 8.99765147e-09 7.30783716e-09 5.54323751e-09 3.72380047e-09
1.86990862e-09 -1.95912605e-23 -1.86990862e-09 -3.72380047e-09
-5.54323751e-09 -7.30783716e-09 -8.99765147e-09]
[ 8.90307127e-09 7.23234206e-09 5.48776061e-09 3.68804375e-09
1.85262574e-09 1.64942482e-23 -1.85262574e-09 -3.68804375e-09
-5.48776061e-09 -7.23234206e-09 -8.90307127e-09]
....
[ 8.99765147e-09 7.30783716e-09 5.54323751e-09 3.72380047e-09
1.86990862e-09 -1.95912605e-23 -1.86990862e-09 -3.72380047e-09
-5.54323751e-09 -7.30783716e-09 -8.99765147e-09]]
ykick : [[ 1.54994519e-08 1.56885817e-08 1.58333555e-08 1.59330783e-08
1.59900727e-08 1.60083620e-08 1.59900727e-08 1.59330783e-08
1.58333555e-08 1.56885817e-08 1.54994519e-08]
....
[-1.54994519e-08 -1.56885817e-08 -1.58333555e-08 -1.59330783e-08
-1.59900727e-08 -1.60083620e-08 -1.59900727e-08 -1.59330783e-08
-1.58333555e-08 -1.56885817e-08 -1.54994519e-08]]
xtable : [-0.01 -0.008 -0.006 -0.004 -0.002 0. 0.002 0.004 0.006 0.008
0.01 ]
ytable : [-0.005 -0.004 -0.003 -0.002 -0.001 0. 0.001 0.002 0.003 0.004
0.005]

print(w150_150mm)
InsertionDevice:
FamName : w150_150mm
Length : 1.5
PassMethod : IdTablePass
Filename_in : kickmap_w150_150mm.dat
Energy : 6.0
Nslice : 150
xkick : [[ 1.68056930e-11 1.36364468e-11 1.03404699e-11 6.94795088e-12
3.49035321e-12 2.51420093e-26 -3.49035321e-12 -6.94795088e-12
-1.03404699e-11 -1.36364468e-11 -1.68056930e-11]
....
[ 1.68056930e-11 1.36364468e-11 1.03404699e-11 6.94795088e-12
3.49035321e-12 2.51420093e-26 -3.49035321e-12 -6.94795088e-12
-1.03404699e-11 -1.36364468e-11 -1.68056930e-11]]
ykick : [[ 4.86245762e-11 4.94306764e-11 5.00658101e-11 5.05238903e-11
5.08005148e-11 5.08930201e-11 5.08005148e-11 5.05238903e-11
5.00658101e-11 4.94306764e-11 4.86245762e-11]
...
[-4.86245762e-11 -4.94306764e-11 -5.00658101e-11 -5.05238903e-11
-5.08005148e-11 -5.08930201e-11 -5.08005148e-11 -5.05238903e-11
-5.00658101e-11 -4.94306764e-11 -4.86245762e-11]]
xtable : [-0.01 -0.008 -0.006 -0.004 -0.002 0. 0.002 0.004 0.006 0.008
0.01 ]
ytable : [-0.005 -0.004 -0.003 -0.002 -0.001 0. 0.001 0.002 0.003 0.004
0.005]

Thank you for your feedback on these tests!

best regards
Simone

[oscarxblanco]

Dear @simoneliuzzo , good to know it is of any use also outside SOLEIL.
The kickmap was a theoretical approximation derived by Pascale ELLEAUME, for the interaction of a parallel electron beam with a potential from magnetic fields (in Eq. 5 the potential phi has units of T² m²)
https://accelconf.web.cern.ch/e92/PDF/EPAC1992_0661.PDF

During the development, however, I did not check inside the IdTablePass method, I only checked that the input to that method was the same in python and MATLAB. Therefore, the 1./Brho**2 is just a factor derived from the attable_id() function in AT MATLAB.

I'm attaching to this message a zip file with : a matlab script where I check the tune variation in matlab, a python script where I check the same in python, the lattice I used, and as you requested the Radia field map file.
example_InsertionDevice_2023feb07.zip

In order to obtain the tune variation, I manually set the InsertionDevice pass method to "DriftPass" for the initial optics, I then set again manually the pass method to "IdTablePass" and recalculate the optics functions. Lastly, I subtract one from the other.

## choose pass method, for a quick test
#newlattice[insertIDindex].PassMethod = 'DriftPass'
#newlattice[insertIDindex].PassMethod = 'IdTablePass'

I guess this procedure is also necessary to get the beta-beat, normalizing where necessary.

I did not, however, set any other parameter that was not used by the IdTablePass. Therefore, the element has no PolynomB or PolynomA properties. I'm not sure how the IdTable pass method is taken into account by the linear optics functions, but, from your plots it seems definitely that something is missing.

In AT MATLAB I saw that the PolynomB and A were set to zero, but, it does produce a perturbation of the optical functions along the ring.

I hope this helps, looking forward to your feedback.

Best regards,
o

[lfarv]

Hello @oscarxblanco and @simoneliuzzo : could you join you kickmaps to your posts, please ?

Concerning energy scaling, the maps should contain in blocks 1 and 2 the ID intrinsic effect, which scale with 1/γ² (so the 1/brho^2 make sense), and in blocks 4 and 5 optional ID errors which scale with 1/γ.

[oscarxblanco]

@lfarv
u20_g45mm_kicks_2022nov10.txt

The errors table (block 4 and 5) are not implemented in element, yet. The factor 1/gamma is written in line 1188 but never used.

Best regards,
o

[simoneliuzzo]

Dear @lfarv as you asked,

I attach the kickmaps provided by Gael to this message.
Notice that the positions are in mm and not in m.

thank you for your help

best regards
Simone

W150.zip

[simoneliuzzo]

Dear @oscarxblanco ,

you are correct! my indexing was correct, but the length was not, and so also the entrance of the elements.

I fix the test with the modification you propose. The figure is now much better.

[oscarxblanco]

@simoneliuzzo very nice to hear.
Do you obtain the same with matlab ?

[oscarxblanco]

@lfarv

Nice addition, thanks
Given the size of the code, I think you could move it to a dedicated module elements.py now really becomes very long…

Is there any naming convention I should follow, or a suggested name for the element and new file name ?

Best regards,
o

[lfarv]

Is there any naming convention I should follow, or a suggested name for the element and new file name ?

Could be an IdTable element and a idtable_element.py file, but up to now there is no strict conventions…

[oscarxblanco]

Dear @lfarv , dear @simoneliuzzo

I have just separated the element in a new file and changed the element name (it is not backwards compatible, I didn'it see the need as this is not released yet).
The new element name is now IdTable as suggested by @lfarv .

u20 = at.IdTable('u20',10,"u20_g45mm_kicks_2022nov10.txt", 2.75)

I profitted to correct the mistakes in my code wrt pep8 code standards, pointed out by @lfarv (thank you again for the tools suggestions). As it is the first time I check for pep8 I only modified the points highlighted by the tool.

Let me know if these changes are OK.

Best regards,
o

[oscarxblanco]

Dear all,
I was thinking that an _activate method to set the pass method to IdTablePass, and a _deactivate method to set the pass method to DriftPass would be useful due to the change in indexing, and elements' length around the ID when inserting the new element.

Let me know if you consider this coherent with the rest of pyat.

Best regards,
o

[simoneliuzzo]

Dear @oscarxblanco and @lfarv ,

a few suggestion for the name, if we are still in time to change:

InsertionDeviceIntegrtedKickMap
InsertionDeviceKickMap
InsertionDeviceKickTable
InsertionDeviceRadiaKickMap

InsertionDeviceIntegrtedKicksMap
InsertionDeviceKicksMap
InsertionDeviceKicksTable
InsertionDeviceRadiaKicksMap

InsertionDeviceIntegrted2DKicksMap
InsertionDevice2DKicksMap
InsertionDevice2DKicksTable
InsertionDeviceRadia2DKicksMap

[oscarxblanco]

@simoneliuzzo
There is no problem in changing it, I will wait a little bit to see if there is any more suggestions.

[simoneliuzzo]

@simoneliuzzo There is no problem in changing it, I will wait a little bit to see if there is any more suggestions.

Thanks, I actually liked much more the original InsertionDevice compared to the present IdTable. An InsertionDevice element may have a pass method IdTable or ID_mega_map or non_linear_id_3D_symplectic_exact_rad_pass, etc...

[lfarv]

@simoneliuzzo, @oscarxblanco:
I agree that a generic InsertionDevice could work, and even the existing Wiggler element. But the required attributes of this new element are strongly related to the presence of a map, whatever it is. There can still be various ways to use a map in different passmethods, but I'd suggest that "Map" appears somewhere in the element name: the InsertionDeviceKickMap that you proposed, for instance?

[oscarxblanco]

Dear @simoneliuzzo , @lfarv ,

I would then opt for the name InsertionDeviceIntegratedKickMaps
as it explicitely refers to an insertion Device modelled by an integrated set of kick maps.

I will implement the change in name, and push it later today. Let me know any other comment/suggestion.
Best regards,
o

[lfarv]

@oscarxblanco: by being so precise in the element name, you forbid any other use of the element with another kind of map and another passmethod… It's rather restrictive.

[oscarxblanco]

@lfarv You are right, and while making changes I saw that my proposal was an awefully long name.
So, I took your suggestions and name it InsertionDeviceKickMap

As an example,

u20 = at.InsertionDeviceKickMap('u20',10,"u20_g45mm_kicks_2022nov10.txt", 2.75)

[simoneliuzzo]

I have used this pull-request branch and it is working well in python

[oscarxblanco]

Dear all,

the element is working well as it is, but, there are functions that will try to split the elements into several pieces in a way that is incompatible with an integrated table. There are two ways out of this problem, as far as I see: either we never split elements with pass method 'IdTablePass' which would allow to avoid assumptions (I.e. if necessary for plots or details in the optics functions one could insert shorter integrated maps); or, to be more flexible, add an scale factor proportional to the length so that the integrated kick is also splitted.

The second option would require more work and care because one could either implement the split as a method of the element, and keep track of the fraction it represents so one should check in every split if there exists already a splitting method, or try to deal with the split in a generic function and loose the connection to the original input file (as the filename, file content, and table values do not reflect the content in memory).

Let me know your opinions.

Best regards,
o

[lfarv]

@oscarxblanco : obviously your 1st option (never split an element with IdTablePass) must be implemented. This should be solved in another pull request, and for me it should not prevent merging this one (anyway the bug is there in Matlab).

[lfarv]

@oscarxblanco: I think this branch can be merged, if you are happy with it.

[oscarxblanco]

@lfarv yes, I think this is ok for a first version. I will leave here a message to recapitulate.

This branch implements a new element, called InsertionDeviceKickMap, in pyat. It uses the pass method IdTablePass for tracking, and requires : a name, the number of slices, an input file with a magnetic field map from Radia, and the beam energy in GeV. For example,

u20 = at.InsertionDeviceKickMap('u20',10,"u20_g45mm_kicks_2022nov10.txt", 2.75)

The InsertionDevice element has been checked against AT MATLAB results for a particular magnetic field map, energy and ring, giving the same tune variation and frequency map.

In addition, the function lattice.insert is incompatible with the InsertionDevice parameters because insert does not copy the Energy property of the element. As a workaround, I have created a dummy element in the lattice that is overwritten by the insertion device.

# create a new lattice where we will insert the InsertionDevice
newlattice = lattice.deepcopy()
### 2023feb05: insert function does not work, it removes the Energy property
### newlattice.insert(insertIDindex+1, u20)
### work around because at.lattice.insert does not copy the Energy property
dummyelem = newlattice[insertIDindex].deepcopy()
newlattice.insert(insertIDindex, dummyelem)
newlattice[insertIDindex] = u20.deepcopy()
# substract half length on each side
newlattice[insertIDindex - 1].Length = newlattice[insertIDindex - 1].Length - u20.Length/2;
newlattice[insertIDindex + 1].Length = newlattice[insertIDindex + 1].Length - u20.Length/2;

In order to obtain the tune variation, I manually set the InsertionDevice pass method to "DriftPass" for the initial optics, I then set again manually the pass method to "IdTablePass" and recalculate the optics functions. Lastly, I subtract one from the other.

## choose pass method, for a quick test
#newlattice[insertIDindex].PassMethod = 'DriftPass'
#newlattice[insertIDindex].PassMethod = 'IdTablePass'

I guess this procedure is also necessary to get the beta-beat, normalizing where necessary.
One could also use the element methods : set_DriftPass, set_IdTablePass, and get_PassMethod, if more convinient. For example,

# newlattice[insertIDindex].set_DriftPass()
newlattice[insertIDindex].set_IdTablePass()
print(newlattice[insertIDindex].get_PassMethod())

NOT INCLUDED :

1. The element has no PolynomB or PolynomA properties
2. The input field map file only contains blocks 1 and 2, The errors table (block 4 and 5) are not implemented in this element, yet. The factor 1/gamma is written in line 1188 but never used.
3. The element does not have yet a function to be read or written, from or to, .mat or .m files.
4. The element cannot be split