Simple ring model

[lcarver]

This is presented as a solution to issue #606 . This PR was also a big motivation in the creation of issue #642 which has not yet been solved.

Here I have introduced a way to create a very simple ring based on a dictionary of global machine parameters. The produced lattice has the following elements.

1. An RF Cavity element based on the provided parameters
2. An M66 matrix that provides tbt tracking for the given Qx and Qy and the path lengthening in z.
3. A detuning element that covers linear chromaticity and amplitude detuning
4. An EnergyLoss element that applies the U0.
5. A new passmethod 'SimpleQuantDiffPass' which applies random noise to xp, yp, and dp with a given sigma that can recreate the desired equilibrium emittances.

I added the new passmethod, added an element to elements.py, and added a function to fastring.py. When we have a solution to #642, we can simply add the damping partition numbers as a parameter in the function 'gen_simple_ring' and pass the damping partition numbers to the EnergyLoss element.

Below is a copy paste of a test file that I have written that shows good agreement for tunes, damping times, equilibrium emittances and z center of mass.

OUTDATED. SEE NEWFILE BELOW

[lcarver]

If it is not possible to easily add damping partition numbers to the EnergyLoss element. I can change a bit how the effects are introduced.

1. I remove the energyloss element.
2. I add the U0 as a misalignment error.
3. I provide damping 'artificially' with user specified damping times and the damping is introduced either in SimpleQuantDiffPass or by adding a damping term to the Mat66.

[swhite2401]

@lcarver, before going into the details of the implementation I would suggest to integrate this directly in the fastring function instead of adding a new one: make the new function private and call it in fastring in case the ring is not provided.

[swhite2401]

Similarly, the quantum diffusion element could be integrate in radiation.gen_quantdiffelem when the ring is not provided

[lfarv]

I fully agree with @swhite2401 on one point: the output should be similar to the one of fastring: exactly the same list of elements with the same functionality. I don't see why you need a new passmethod if a fastring can represent everything.

On the other hand, I don't see a good reason to group both in a same function: the inputs are completely different, and I don't think there is anything in common in the implementation. fastring and simplering look nice for me.

[swhite2401]

A passmethod to model damping and diffusion with only given equilibrium emittance and damping times is a useful addition. It allows to generate a fastring from a parameter table when the full AT lattice file is not available which is the whole point of this PR.

The 3 steps that @lcarver proposed are fine and equivalent to what is done in fastring.

I disagree that the structure/functionaltity should be the same as fastring, I can see some advantage in having all synchrotron radiation effects in a single element and not mixed with the M66 as is done in fastring. This allows to turn on/off radiation very easily and prevents from having to generate 2 separate ring with and w.o. radiations as we presently do with the fast ring.

[lfarv]

I have no objection to split betatron motion and damping in two different elements, but:

• then one should do the same in fastring: your argument about turning on/off is equally valid for fastring. It's trivial to split the global M66 matrix,
• damping (excluding diffusion) does not require a new pass method: M66 is ok for that, and generating the matrix from the given damping times is easy.
• diffusion and energy loss can be treated exactly as in fastring.

But still none of the proposed splits allows emittance computation with ohmienvelope… Too bad.

[swhite2401]

I prefer what is proposed by @lcarver, is it much simpler and more efficient. Having all radiation effects in a single element is a very nice and convenient feature and even if it is not used for this simple fastring this new passmethd is a nice addition that should be integrated in AT.

[lfarv]

Why not, but @lcarver proposes to include the damping in the (single) 6x6 matrix, while I thought you preferred having it separated from betatron motion.

Also, the pass method included in this PR does not deal neither with damping nor with energy loss: it only does the quantum excitation.

Treating damping and energy loss is the "classical" way has some advantages: if you disable it, the beams damps to zero but the closed orbit and linear optics still give the good results: you get the right synchronous phase and damping times.

Anyway, I don't care to much about how the transfer is split, I just think that once the optimum is chosen, there is no reason why it should be different between fastring and simplering…

[carmignani]

I think it would be nice to have the damping separated from the M66 in the fast ring, so that we could simply remove the diffusion element and the damping element and have the fast ring without radiation. The simplest way is to have an additional linear matrix doing only the damping, with the disadvantage that the tracking would be slower. Having a single synchrotron radiation element doing both damping and diffusion wouldn't allow to do simulations with only damping, that is what we often do in single particle simulations. The diffusion matrix for the simple ring could be built from the moments matrix and the M66 matrix (with damping) without the need of an additional passmethod.

[lcarver]

I wrote a comment here but I think I had a wrong interpretation of some things.

I will implement a correct (without EnergyLoss element) simple ring model, that includes a SimpleQuantumDiffusion element that combines both radiation damping and equilibrium emittance. By lumping into one element, this will for sure be faster. I will include a more detailed comment when it is ready for review.

[lcarver]

The SimpleQuantDiffPass has been modified. It now requires U0 (in addition to previous arguments). The passmethod now implements radiation damping (if tau>0) and the simple quantum diffusion (if emit>0) and it also removes U0/ring->Energy (same as applying a misalignment) to model the energy loss.

If the user asks for emit=0, then no quantum diffusion is applied. If the user asks for emit=0 and tau=0, then radiation damping and quantum diffusion is not applied

The energyloss element has been removed as U0 is covered in SimpleQuantDiffPass with user defined damping times. The file below is the updated one from the first comment (which I have now removed). It shows good damping times and emittances, correct tunes and beam centered on zero.

For me this solves the problem of generating a ring from a table. This is the simplest implementation and is a very straightforward way of generating a lattice that can be used for collective effects studies in the fewest elements. Ready for review.

import at
from at.physics.fastring import gen_simple_ring
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import hilbert
from at.constants import clight

dic = {
       'alpha_x': 1e-3,
       'alpha_y': 1e-4,
       'beta_x': 20,
       'beta_y': 10,
       'Qx': 74.1235,
       'Qy': 28.3,
       'Qpx': 0,
       'Qpy': 0,
       'A1': 0.,
       'A2': 0.,
       'A3': 0.,
       'circumference': 843.977,
       'alpha': 8.5e-5,
       'U0': 2.5e6,
       'Vrf': 6e6,
       'energy': 6e9,
       'harmonic_number': 992,
       'emit_x': 100e-12,
       'emit_y': 10e-12,
       'sigma_dp': 9.5e-4,
       'tau': [3000, 5000, 1000]
}

ring = gen_simple_ring(dic)

r_in = np.zeros((6,1))
r_in[0] = 1e-1
r_in[2] = 1e-1
r_in[4] = 1e-3

# Check radiation damping
nt = 10000
r_out = at.lattice_pass(ring, r_in, nturns=nt)

xdat = r_out[0,0,0,:]
fitx = np.polyfit(np.arange(2500), np.log(np.abs(hilbert(xdat)))[500:3000], 1)
print('tau_x', -1/fitx[0], dic['tau'][0])

ydat = r_out[2,0,0,:]
fity = np.polyfit(np.arange(2500), np.log(np.abs(hilbert(ydat)))[500:3000], 1)
print('tau_y', -1/fity[0], dic['tau'][1])

#spikes at start and finish are artefacts of hilbert transform
plt.suptitle('envelope from hilbert transform')
plt.semilogy(np.abs(hilbert(xdat)),'r', label='x')
plt.semilogy(np.abs(hilbert(ydat)),'b', label='y')
plt.xlabel('turn')
plt.ylabel('env [m]')
plt.legend()
plt.show()

r_in = np.zeros((6,1))
r_in[0] = 0.
r_in[2] = 0.
r_in[4] = 1e-3
r_out = at.lattice_pass(ring, r_in, nturns=nt)
r_in[4] += 5e-2
r_out = at.lattice_pass(ring, r_in, nturns=nt)
dpdat = r_out[4,0,0,:]
zdat = r_out[5,0,0,:]

fitz = np.polyfit(np.arange(500), np.log(np.abs(hilbert(dpdat)))[250:750], 1)
print('tau_z', -1/fitz[0], dic['tau'][2])

plt.suptitle('abs(hilbert(dp))')
plt.semilogy(np.abs(hilbert(dpdat)))
plt.show()

# Check emittances
sigma_matrix = at.sigma_matrix(betax=dic['beta_x'], betay=dic['beta_y'], emitx=dic['emit_x'], emity=dic['emit_y'], espread=dic['sigma_dp'], blength=3e-4, alphax=dic['alpha_x'], alphay=dic['alpha_y'])
nparts = 10000

r_in = at.beam(nparts, sigma_matrix)

emitx = []
emity = []
sigma_e = []
blength = []
for i in np.arange(20000):
    _ = at.lattice_pass(ring, r_in, nturns=1)

    track_emitx = np.sqrt(np.mean(r_in[0,:]**2)*np.mean(r_in[1,:]**2) - np.mean(r_in[0,:]*r_in[1,:])**2)
    track_emity = np.sqrt(np.mean(r_in[2,:]**2)*np.mean(r_in[3,:]**2) - np.mean(r_in[2,:]*r_in[3,:])**2)
    emitx.append(track_emitx)
    emity.append(track_emity)
    sigma_e.append(np.std(r_in[4,:]))
    blength.append(np.std(r_in[5,:]))

print('emitx', track_emitx, dic['emit_x'])
print('emity', track_emity, dic['emit_y'])
print('sigma_e', np.std(r_in[4,:]), dic['sigma_dp'])
print('blength', np.std(r_in[5,:]))

# Check tunes
freq = np.fft.rfftfreq(nt)
ffreq = np.fft.rfftfreq(nt,d=dic['circumference']/clight)
plt.semilogy(freq, np.abs(np.fft.rfft(xdat)),'r', label='x')
plt.axvline(dic['Qx']%1, color='r')

plt.semilogy(freq, np.abs(np.fft.rfft(ydat)), 'b', label='y')
plt.axvline(dic['Qy']%1, color='b')
plt.legend()
plt.xlabel('tune')
plt.ylabel('abs(fft)')
plt.show()

#Expected fs is about 1300 Hz.
plt.semilogy(ffreq, np.abs(np.fft.rfft(dpdat)), 'b')
plt.xlim([0,3000])
plt.xlabel('freq [Hz]')
plt.ylabel('abs(fft(dp))')
plt.show()

ring_noqd = ring.deepcopy()
u0_test = ring_noqd.get_energy_loss(at.energy_loss.ELossMethod.TRACKING)
print(u0_test)

# Check TimeLag of cavity
r_in = np.zeros((6,1))
r_out = at.lattice_pass(ring, r_in, nturns=100000)
print('zc', np.mean(r_out[5,0,0,:]))

[lcarver]

Hello, All of your modifications were taken into account. Default passmethod of SimpleQuantDiff is SimpleQuantDiffPass. When applying rad off it switches to IdentityPass.

Rather than implement the optional / default values in the C passmethod, I prefer to set them all on the python side. this allows me to make some assertions about what cases are allowed or not (tau with emit=0 is ok, emit with tau=0 is not).

Now, simple_ring requires 7 positional arguments: energy, circumference, harmonic_number, Qx, Qy, Vrf, alpha. The rest of the parameters are keyword arguments that have default values.

New parsing of the function looks like this:

dic = {
       'alpha_x': 1e-3,
       'alpha_y': 1e-4,
       'beta_x': 20,
       'beta_y': 10,
       'Qx': 74.1235,
       'Qy': 28.3,
       'Qpx': 0,
       'Qpy': 0,
       'A1': 0.,
       'A2': 0.,
       'A3': 0.,
       'circumference': 843.977,
       'alpha': 8.5e-5,
       'U0': 2.5e6,
       'Vrf': 6e6,
       'energy': 6e9,
       'harmonic_number': 992,
       'emit_x': 100e-12,
       'emit_y': 10e-12,
       'sigma_dp': 9.5e-4,
       'tau_x': 3000,
       'tau_y': 5000,
       'tau_z': 1000
}

ring = simple_ring(dic['energy'], dic['circumference'], dic['harmonic_number'], dic['Qx'], dic['Qy'], dic['Vrf'], dic['alpha'], tau_x=dic['tau_x'], emit_x=dic['emit_x'], tau_y=dic['tau_y'], emit_y=dic['emit_y'], tau_z=dic['tau_z'], sigma_dp=dic['sigma_dp'], U0=dic['U0'], A1=dic['A1'], A2=dic['A2'], A3=dic['A3'], Qpx=dic['Qpx'], Qpy=dic['Qpy'])

but similarly you can just use

ring = simple_ring(6e9, 844, 992, 0.1, 0.1, 6e6, 8.5e-5)

[lcarver]

ready for review

[lcarver]

all implemented except the last comment