From 2a3d5e2e04a45aaa46ee358ee3d006c15cc51cc8 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 9 Dec 2024 10:55:06 +0100
Subject: [PATCH 01/48] changed indexing in Ring .delta_E to work for multiple
 sections

---
 blond/input_parameters/ring.py | 15 ++++++++++++++-
 1 file changed, 14 insertions(+), 1 deletion(-)

diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 3040f51a..57fb16a5 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -277,12 +277,25 @@ class Ring:
         self.energy: NDArray = np.sqrt(self.momentum ** 2 + self.particle.mass ** 2)
         self.kin_energy: NDArray = (np.sqrt(self.momentum ** 2 + self.particle.mass ** 2) -
                                        self.particle.mass)
-        self.delta_E: NDArray = np.diff(self.energy, axis=1)
         self.t_rev: NDArray = np.dot(self.ring_length, 1 / (self.beta * c))
         self.cycle_time: NDArray = np.cumsum(self.t_rev)  # Always starts with zero
         self.f_rev: NDArray = 1 / self.t_rev
         self.omega_rev: NDArray = 2 * np.pi * self.f_rev
 
+
+                
+        if self.n_sections == 1:
+            self.delta_E = np.diff(self.energy, axis=1)
+        else:
+            self.delta_E = np.reshape(np.zeros(n_turns*n_sections),(n_sections,n_turns))
+            for i in range(n_turns):  
+                for j in range(n_sections):
+                    if j == 0:
+                        self.delta_E[j,i] = self.energy[j,i+1]-self.energy[-1,i]
+                    else:
+                        self.delta_E[j,i] = self.energy[j,i+1]-self.energy[j-1,i+1]
+
+
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
             interp_time = self.cycle_time
-- 
GitLab


From d5ef5a026d1445aa5b255b69e974f7edd0e175d0 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Thu, 12 Dec 2024 15:46:25 +0100
Subject: [PATCH 02/48] changed variable inheritance for
 InducedVoltageResonator

---
 blond/impedances/impedance.py | 36 ++++++++++++++++++++++++-----------
 1 file changed, 25 insertions(+), 11 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index e1a29d16..de0e5897 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -340,10 +340,11 @@ class _InducedVoltage:
         if self.multi_turn_wake:
             # Number of points of the memory array for multi-turn wake
             self.n_mtw_memory = self.n_induced_voltage
-
+            print('mtw initialised')
             self.front_wake_buffer = 0
 
             if self.mtw_mode == 'freq':
+                print('in freq domain')
                 # In frequency domain, an extra buffer for a revolution turn is
                 # needed due to the circular time shift in frequency domain
                 self.buffer_size = np.ceil(np.max(self.rf_params.t_rev)
@@ -365,6 +366,7 @@ class _InducedVoltage:
                 # Selecting time-shift method
                 self.shift_trev = self.shift_trev_freq
             else:
+                print('in time domain')
                 # Selecting time-shift method
                 self.shift_trev = self.shift_trev_time
                 # Time array
@@ -430,7 +432,7 @@ class _InducedVoltage:
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
 
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
-
+        print('multi turn wake is being calculated')
     def shift_trev_freq(self):
         """
         Method to shift the induced voltage by a revolution period in the
@@ -946,10 +948,20 @@ class InducedVoltageResonator(_InducedVoltage):
         Computed induced voltage [V]
     """
 
-    @handle_legacy_kwargs
-    def __init__(self, beam: Beam, profile: Profile, resonators: Resonators,
-                 timeArray: Optional[NDArray] = None):
-
+        
+    @handle_legacy_kwargs 
+    def __init__(self, beam: Beam,
+                 profile: Profile,
+                 resonators: list[_ImpedanceObject],
+                 frequency_resolution: Optional[float] = None,
+                 wake_length: Optional[float] = None,
+                 multi_turn_wake: bool = False,
+                 timeArray: Optional[NDArray] = None,
+                 rf_station: Optional[RFStation] = None,
+                 mtw_mode: Optional[MtwModeTypes] = None,
+                 use_regular_fft: bool = True) -> None:
+        
+        
         # Test if one or more quality factors is smaller than 0.5.
         if sum(resonators.Q < 0.5) > 0:
             # ResonatorError
@@ -1000,11 +1012,13 @@ class InducedVoltageResonator(_InducedVoltage):
         self._deltaT = np.zeros(
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
-        # Call the __init__ method of the parent class [calls process()]
-        _InducedVoltage.__init__(self, beam, profile, wake_length=None,
-                                 frequency_resolution=None,
-                                 multi_turn_wake=False, rf_station=None,
-                                 mtw_mode=None)
+         # Call the __init__ method of the parent class [calls process()]
+        _InducedVoltage.__init__(self, beam, profile,
+                                 frequency_resolution=frequency_resolution,
+                                 wake_length=wake_length,
+                                 multi_turn_wake=multi_turn_wake,
+                                 rf_station=rf_station, mtw_mode=mtw_mode,
+                                 use_regular_fft=use_regular_fft)
 
     def process(self):
         r"""
-- 
GitLab


From e384976f5968e331c67e865bda6bee16c7e60b7d Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Mon, 16 Dec 2024 11:29:35 +0100
Subject: [PATCH 03/48] Wrote test case stubs

---
 unittests/impedances/test_impedance.py | 74 ++++++++++++++++++++++++--
 1 file changed, 70 insertions(+), 4 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index ce73d653..541aefa5 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -18,14 +18,13 @@ import unittest
 import numpy as np
 
 from blond.beam.profile import CutOptions, Profile
-from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime
+from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime, InducedVoltageResonator
 from blond.impedances.impedance_sources import Resonators
 
 
 class TestInducedVoltageFreq(unittest.TestCase):
 
     def setUp(self):
-
         # Nyquist frequency 1.6 GHz; frequency spacing 50 MHz
         self.profile = Profile(None,
                                cut_options=CutOptions(cut_left=0, cut_right=5e-9, n_slices=16))
@@ -86,7 +85,6 @@ class TestInducedVoltageFreq(unittest.TestCase):
 class TestInducedVoltageTime(unittest.TestCase):
 
     def setUp(self):
-
         # Nyquist frequency 1.6 GHz; frequency spacing 50 MHz
         self.profile = Profile(None,
                                cut_options=CutOptions(cut_left=0, cut_right=5e-9, n_slices=16))
@@ -121,6 +119,74 @@ class TestInducedVoltageTime(unittest.TestCase):
         np.testing.assert_allclose(test_object.wake_length_input, 11e-9)
 
 
-if __name__ == '__main__':
+class TestInducedVoltageResonator(unittest.TestCase):
 
+    def setUp(self):
+        import os
+        import numpy as np
+        from blond.beam.beam import Beam, Proton
+        from blond.beam.distributions import bigaussian
+        from blond.beam.profile import CutOptions, FitOptions, Profile
+        from blond.impedances.impedance_sources import Resonators
+        from blond.input_parameters.rf_parameters import RFStation
+        from blond.input_parameters.ring import Ring
+        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+        gamma_transition = 1 / np.sqrt(0.00192)  # [1]
+        ring = Ring(
+            ring_length=6911.56,
+            alpha_0=1 / gamma_transition ** 2,
+            synchronous_data=25.92e9,
+            particle=Proton(),
+            n_turns=2,
+        )
+        rf_station = RFStation(
+            ring=ring,
+            harmonic=[4620],
+            voltage=[0.9e6],
+            phi_rf_d=[0.0],
+            n_rf=1,
+        )
+        beam = Beam(
+            ring=ring,
+            n_macroparticles=1001,
+            intensity=int(1e10),
+        )
+        bigaussian(
+            ring=ring,
+            rf_station=rf_station,
+            beam=beam,
+            sigma_dt=2e-9 / 4, seed=1,
+        )
+        self.beam = beam
+        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
+                                 rf_station=rf_station, cuts_unit='rad')
+        self.profile = Profile(beam, cut_options,
+                               FitOptions(fit_option='gaussian'))
+        table = np.loadtxt(this_directory + '/EX_05_new_HQ_table.dat', comments='!')
+        R_shunt = table[:, 2] * 10 ** 6
+        f_res = table[:, 0] * 10 ** 9
+        Q_factor = table[:, 1]
+        self.resonator = Resonators(R_shunt, f_res, Q_factor)
+
+    def test_init(self):
+        # TODO Improve testcases
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
+
+    def test_init_mtw(self):
+        # TODO Improve testcases
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      mtw_mode=True
+                                      )
+
+    def test_init_mtw2(self):
+        # TODO Improve testcases
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      mtw_mode=True
+                                      )
+        ivr.process()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+
+if __name__ == '__main__':
     unittest.main()
-- 
GitLab


From e818d9e7169d18b67544b41a208f3c9ad91f510a Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Mon, 16 Dec 2024 11:48:54 +0100
Subject: [PATCH 04/48] mtw working testcase

---
 blond/impedances/impedance.py          | 53 ++++++++++++++------------
 unittests/impedances/test_impedance.py | 12 +++---
 2 files changed, 36 insertions(+), 29 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index de0e5897..3df43d09 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -433,6 +433,7 @@ class _InducedVoltage:
 
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
         print('multi turn wake is being calculated')
+
     def shift_trev_freq(self):
         """
         Method to shift the induced voltage by a revolution period in the
@@ -534,12 +535,15 @@ class InducedVoltageTime(_InducedVoltage):
         ####################################
 
         # Call the __init__ method of the parent class [calls process()]
-        _InducedVoltage.__init__(self, beam, profile,
-                                 frequency_resolution=None,
-                                 wake_length=wake_length,
-                                 multi_turn_wake=multi_turn_wake,
-                                 rf_station=rf_station, mtw_mode=mtw_mode,
-                                 use_regular_fft=use_regular_fft)
+        super().__init__(beam=beam,
+                         profile=profile,
+                         frequency_resolution=None,
+                         wake_length=wake_length,
+                         multi_turn_wake=multi_turn_wake,
+                         rf_station=rf_station,
+                         mtw_mode=mtw_mode,
+                         use_regular_fft=use_regular_fft,
+                         )
 
     def process(self) -> None:
         """
@@ -714,11 +718,11 @@ class InducedVoltageFreq(_InducedVoltage):
         ###############
 
         # Call the __init__ method of the parent class
-        _InducedVoltage.__init__(self, beam, profile, wake_length=None,
-                                 frequency_resolution=frequency_resolution,
-                                 multi_turn_wake=multi_turn_wake,
-                                 rf_station=rf_station, mtw_mode=mtw_mode,
-                                 use_regular_fft=use_regular_fft)
+        super().__init__(beam, profile, wake_length=None,
+                         frequency_resolution=frequency_resolution,
+                         multi_turn_wake=multi_turn_wake,
+                         rf_station=rf_station, mtw_mode=mtw_mode,
+                         use_regular_fft=use_regular_fft)
 
     def process(self) -> None:
         """
@@ -850,7 +854,7 @@ class InductiveImpedance(_InducedVoltage):
         self.deriv_mode = deriv_mode
 
         # Call the __init__ method of the parent class
-        _InducedVoltage.__init__(self, beam, profile, rf_station=rf_station)
+        super().__init__(beam, profile, rf_station=rf_station)
 
     def induced_voltage_1turn(self, beam_spectrum_dict={}):
         """
@@ -948,11 +952,10 @@ class InducedVoltageResonator(_InducedVoltage):
         Computed induced voltage [V]
     """
 
-        
-    @handle_legacy_kwargs 
+    @handle_legacy_kwargs
     def __init__(self, beam: Beam,
                  profile: Profile,
-                 resonators: list[_ImpedanceObject],
+                 resonators: Resonators,
                  frequency_resolution: Optional[float] = None,
                  wake_length: Optional[float] = None,
                  multi_turn_wake: bool = False,
@@ -960,8 +963,7 @@ class InducedVoltageResonator(_InducedVoltage):
                  rf_station: Optional[RFStation] = None,
                  mtw_mode: Optional[MtwModeTypes] = None,
                  use_regular_fft: bool = True) -> None:
-        
-        
+
         # Test if one or more quality factors is smaller than 0.5.
         if sum(resonators.Q < 0.5) > 0:
             # ResonatorError
@@ -1012,13 +1014,13 @@ class InducedVoltageResonator(_InducedVoltage):
         self._deltaT = np.zeros(
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
-         # Call the __init__ method of the parent class [calls process()]
-        _InducedVoltage.__init__(self, beam, profile,
-                                 frequency_resolution=frequency_resolution,
-                                 wake_length=wake_length,
-                                 multi_turn_wake=multi_turn_wake,
-                                 rf_station=rf_station, mtw_mode=mtw_mode,
-                                 use_regular_fft=use_regular_fft)
+        # Call the __init__ method of the parent class [calls process()]
+        super().__init__(beam, profile,
+                         frequency_resolution=frequency_resolution,
+                         wake_length=wake_length,
+                         multi_turn_wake=multi_turn_wake,
+                         rf_station=rf_station, mtw_mode=mtw_mode,
+                         use_regular_fft=use_regular_fft)
 
     def process(self):
         r"""
@@ -1034,6 +1036,9 @@ class InducedVoltageResonator(_InducedVoltage):
         self._deltaT = np.zeros(
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
+    def on_profile_change(self):
+        self.process(self)
+
     def induced_voltage_1turn(self, beam_spectrum_dict={}):
         r"""
         Method to calculate the induced voltage through linearily
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 541aefa5..13f87e0b 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -139,7 +139,7 @@ class TestInducedVoltageResonator(unittest.TestCase):
             particle=Proton(),
             n_turns=2,
         )
-        rf_station = RFStation(
+        self.rf_station = RFStation(
             ring=ring,
             harmonic=[4620],
             voltage=[0.9e6],
@@ -153,13 +153,13 @@ class TestInducedVoltageResonator(unittest.TestCase):
         )
         bigaussian(
             ring=ring,
-            rf_station=rf_station,
+            rf_station=self.rf_station,
             beam=beam,
             sigma_dt=2e-9 / 4, seed=1,
         )
         self.beam = beam
         cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
-                                 rf_station=rf_station, cuts_unit='rad')
+                                 rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
         table = np.loadtxt(this_directory + '/EX_05_new_HQ_table.dat', comments='!')
@@ -175,15 +175,17 @@ class TestInducedVoltageResonator(unittest.TestCase):
     def test_init_mtw(self):
         # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      mtw_mode=True
+                                      multi_turn_wake=True
                                       )
 
     def test_init_mtw2(self):
         # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      mtw_mode=True
+                                      rf_station=self.rf_station,
+                                      multi_turn_wake=True
                                       )
         ivr.process()
+        ivr.induced_voltage_mtw()
         my_array = ivr.induced_voltage
         self.assertTrue(np.any(my_array != 0.0))
 
-- 
GitLab


From 857c2a19a61c60d8c30cb2ad0766881c4db382d6 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 12:42:56 +0100
Subject: [PATCH 05/48] editing tests

---
 unittests/impedances/test_impedance.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 13f87e0b..b79562f7 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -162,7 +162,7 @@ class TestInducedVoltageResonator(unittest.TestCase):
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        table = np.loadtxt(this_directory + '/EX_05_new_HQ_table.dat', comments='!')
+        table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
         R_shunt = table[:, 2] * 10 ** 6
         f_res = table[:, 0] * 10 ** 9
         Q_factor = table[:, 1]
-- 
GitLab


From f798623e1c88bbc0d9a2a43941459b6e4b10d11c Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 12:43:58 +0100
Subject: [PATCH 06/48] removed print

---
 blond/impedances/impedance.py | 3 ---
 1 file changed, 3 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 3df43d09..65ad49a9 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -340,11 +340,9 @@ class _InducedVoltage:
         if self.multi_turn_wake:
             # Number of points of the memory array for multi-turn wake
             self.n_mtw_memory = self.n_induced_voltage
-            print('mtw initialised')
             self.front_wake_buffer = 0
 
             if self.mtw_mode == 'freq':
-                print('in freq domain')
                 # In frequency domain, an extra buffer for a revolution turn is
                 # needed due to the circular time shift in frequency domain
                 self.buffer_size = np.ceil(np.max(self.rf_params.t_rev)
@@ -366,7 +364,6 @@ class _InducedVoltage:
                 # Selecting time-shift method
                 self.shift_trev = self.shift_trev_freq
             else:
-                print('in time domain')
                 # Selecting time-shift method
                 self.shift_trev = self.shift_trev_time
                 # Time array
-- 
GitLab


From 534acaeac93340d77e5115215178c8f5c2b0022d Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 13:37:36 +0100
Subject: [PATCH 07/48] changed impedance test cases for multi rf
 InducedVoltageResonator

---
 unittests/impedances/test_impedance.py | 149 +++++++++++++++++++++----
 1 file changed, 125 insertions(+), 24 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index b79562f7..88578c77 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -119,7 +119,52 @@ class TestInducedVoltageTime(unittest.TestCase):
         np.testing.assert_allclose(test_object.wake_length_input, 11e-9)
 
 
-class TestInducedVoltageResonator(unittest.TestCase):
+class  BaseTestInducedVoltageResonator(unittest.TestCase):
+
+    def test_init(self):
+        # TODO Improve testcases
+        
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
+
+    def test_init_mtw(self):
+        # TODO Improve testcases
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=True
+                                      )
+        
+    def test_mtw_false_induced_volt(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=False)
+    
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+
+    def test_multi_rf_station(self):
+    
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
+                                      rf_station=self.rf_station, multi_turn_wake=False)
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+        
+
+    def test_mtw_calculation(self):
+        # TODO Improve testcases
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      rf_station=self.rf_station, multi_turn_wake=True)
+        ivr.process()
+        ivr.induced_voltage_mtw()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+
+
+class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.TestCase):
 
     def setUp(self):
         import os
@@ -130,15 +175,23 @@ class TestInducedVoltageResonator(unittest.TestCase):
         from blond.impedances.impedance_sources import Resonators
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
+        import math 
         this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        gamma_transition = 1 / np.sqrt(0.00192)  # [1]
+        
+        gamma_transition = 1 / np.sqrt(0.00192)  
+        mom_compaction = 1 / gamma_transition ** 2
+        circ = 6911.56
+        
         ring = Ring(
-            ring_length=6911.56,
-            alpha_0=1 / gamma_transition ** 2,
+            ring_length=circ,
+            alpha_0=mom_compaction,
             synchronous_data=25.92e9,
             particle=Proton(),
-            n_turns=2,
-        )
+            n_turns=2)
+                        
+        
+        
+
         self.rf_station = RFStation(
             ring=ring,
             harmonic=[4620],
@@ -146,6 +199,7 @@ class TestInducedVoltageResonator(unittest.TestCase):
             phi_rf_d=[0.0],
             n_rf=1,
         )
+        
         beam = Beam(
             ring=ring,
             n_macroparticles=1001,
@@ -168,26 +222,73 @@ class TestInducedVoltageResonator(unittest.TestCase):
         Q_factor = table[:, 1]
         self.resonator = Resonators(R_shunt, f_res, Q_factor)
 
-    def test_init(self):
-        # TODO Improve testcases
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
 
-    def test_init_mtw(self):
-        # TODO Improve testcases
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=True
-                                      )
 
-    def test_init_mtw2(self):
-        # TODO Improve testcases
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station,
-                                      multi_turn_wake=True
-                                      )
-        ivr.process()
-        ivr.induced_voltage_mtw()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
+
+class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.TestCase):
+
+    def setUp(self):
+        import os
+        import numpy as np
+        from blond.beam.beam import Beam, Proton
+        from blond.beam.distributions import bigaussian
+        from blond.beam.profile import CutOptions, FitOptions, Profile
+        from blond.impedances.impedance_sources import Resonators
+        from blond.input_parameters.rf_parameters import RFStation
+        from blond.input_parameters.ring import Ring
+        import math 
+        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+        
+        gamma_transition = 1 / np.sqrt(0.00192)  
+        mom_compaction = 1 / gamma_transition ** 2
+        circ = 6911.56
+        n_rf_stations = 2
+        l_per_section = circ/ n_rf_stations
+
+        ring_per_section = Ring(ring_length=l_per_section, 
+                                alpha_0=mom_compaction,
+                        particle=Proton(), synchronous_data=25.92e9)
+
+        sectionlengths = np.full(n_rf_stations, l_per_section)
+        alpha_c = np.full(n_rf_stations, mom_compaction)
+        n_turns = math.ceil(ring_per_section.n_turns /n_rf_stations)
+        ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c, particle=Proton(), n_turns=n_turns, n_sections=n_rf_stations,
+                    synchronous_data=25.92e9)
+
+
+
+        self.rf_station = RFStation(
+            ring=ring,
+            harmonic=[4620],
+            voltage=[0.9e6],
+            phi_rf_d=[0.0],
+            n_rf=1,
+        )
+        
+        beam = Beam(
+            ring=ring,
+            n_macroparticles=1001,
+            intensity=int(1e10),
+        )
+        bigaussian(
+            ring=ring,
+            rf_station=self.rf_station,
+            beam=beam,
+            sigma_dt=2e-9 / 4, seed=1,
+        )
+        self.beam = beam
+        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
+                                 rf_station=self.rf_station, cuts_unit='rad')
+        self.profile = Profile(beam, cut_options,
+                               FitOptions(fit_option='gaussian'))
+        table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
+        R_shunt = table[:, 2] * 10 ** 6
+        f_res = table[:, 0] * 10 ** 9
+        Q_factor = table[:, 1]
+        self.resonator = Resonators(R_shunt, f_res, Q_factor)
+             
+ 
+
 
 
 if __name__ == '__main__':
-- 
GitLab


From e20f3a913582911160f27a9ec64da26df3eb9d22 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 13:38:23 +0100
Subject: [PATCH 08/48] removed print in impedances.py

---
 blond/impedances/impedance.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 65ad49a9..99deb89d 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -429,7 +429,6 @@ class _InducedVoltage:
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
 
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
-        print('multi turn wake is being calculated')
 
     def shift_trev_freq(self):
         """
-- 
GitLab


From d5119c81b14b7559d538b2a4c25088956f3c7995 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 13:53:22 +0100
Subject: [PATCH 09/48] editing tests

---
 unittests/impedances/test_impedance.py | 34 ++++++++++++++------------
 1 file changed, 19 insertions(+), 15 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 88578c77..ad14cdf5 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -153,7 +153,7 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
         
 
-    def test_mtw_calculation(self):
+    def test_mtw_true_induced_voltage(self):
         # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       rf_station=self.rf_station, multi_turn_wake=True)
@@ -164,6 +164,15 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
 
 
+
+        
+
+                                      
+                        
+
+
+
+
 class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.TestCase):
 
     def setUp(self):
@@ -190,8 +199,6 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
             n_turns=2)
                         
         
-        
-
         self.rf_station = RFStation(
             ring=ring,
             harmonic=[4620],
@@ -217,10 +224,10 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
         table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
-        R_shunt = table[:, 2] * 10 ** 6
-        f_res = table[:, 0] * 10 ** 9
-        Q_factor = table[:, 1]
-        self.resonator = Resonators(R_shunt, f_res, Q_factor)
+        self.R_shunt = table[:, 2] * 10 ** 6
+        self.f_res = table[:, 0] * 10 ** 9
+        self.Q_factor = table[:, 1]
+        self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
 
 
@@ -244,7 +251,6 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
         circ = 6911.56
         n_rf_stations = 2
         l_per_section = circ/ n_rf_stations
-
         ring_per_section = Ring(ring_length=l_per_section, 
                                 alpha_0=mom_compaction,
                         particle=Proton(), synchronous_data=25.92e9)
@@ -252,11 +258,10 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
         sectionlengths = np.full(n_rf_stations, l_per_section)
         alpha_c = np.full(n_rf_stations, mom_compaction)
         n_turns = math.ceil(ring_per_section.n_turns /n_rf_stations)
+
         ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c, particle=Proton(), n_turns=n_turns, n_sections=n_rf_stations,
                     synchronous_data=25.92e9)
 
-
-
         self.rf_station = RFStation(
             ring=ring,
             harmonic=[4620],
@@ -282,14 +287,13 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
         table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
-        R_shunt = table[:, 2] * 10 ** 6
-        f_res = table[:, 0] * 10 ** 9
-        Q_factor = table[:, 1]
-        self.resonator = Resonators(R_shunt, f_res, Q_factor)
+        self.R_shunt = table[:, 2] * 10 ** 6
+        self.f_res = table[:, 0] * 10 ** 9
+        self.Q_factor = table[:, 1]
+        self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
              
  
 
 
-
 if __name__ == '__main__':
     unittest.main()
-- 
GitLab


From 7b75a5cb0f68e3300f64c2c32f17a54fe9d11027 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 14:29:00 +0100
Subject: [PATCH 10/48] remove parameter table from test folder

---
 unittests/impedances/test_impedance.py | 74 ++++++++++++++++++--------
 1 file changed, 53 insertions(+), 21 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index ad14cdf5..60b5c679 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -123,7 +123,7 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
     def test_init(self):
         # TODO Improve testcases
-        
+    
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
 
     def test_init_mtw(self):
@@ -133,6 +133,7 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
                                       )
         
     def test_mtw_false_induced_volt(self):
+        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       multi_turn_wake=False)
     
@@ -143,7 +144,7 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
 
     def test_multi_rf_station(self):
-    
+        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
                                       rf_station=self.rf_station, multi_turn_wake=False)
         ivr.process()
@@ -165,18 +166,9 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
 
 
-        
-
-                                      
-                        
-
-
-
-
 class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.TestCase):
 
     def setUp(self):
-        import os
         import numpy as np
         from blond.beam.beam import Beam, Proton
         from blond.beam.distributions import bigaussian
@@ -185,7 +177,6 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
         import math 
-        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
         
         gamma_transition = 1 / np.sqrt(0.00192)  
         mom_compaction = 1 / gamma_transition ** 2
@@ -223,10 +214,32 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
-        self.R_shunt = table[:, 2] * 10 ** 6
-        self.f_res = table[:, 0] * 10 ** 9
-        self.Q_factor = table[:, 1]
+        self.R_shunt = np.array([3.88e+05, 6.00e+02, 6.98e+04, 8.72e+04, 2.04e+04, 1.09e+05,
+       9.10e+03, 1.80e+03, 7.50e+03, 2.60e+03, 2.80e+04, 2.90e+04,
+       1.20e+05, 7.90e+04, 4.80e+04, 1.45e+04, 5.00e+03, 1.05e+04,
+       1.25e+04, 1.00e+04, 4.30e+04, 4.50e+04, 1.15e+04, 4.50e+03,
+       6.90e+03, 1.54e+04, 1.10e+04, 4.70e+03, 8.50e+03, 1.70e+04,
+       5.00e+03, 6.33e+05, 4.99e+05, 3.64e+05, 1.77e+05, 6.14e+04,
+       2.99e+04, 3.79e+05, 2.43e+05, 1.74e+04, 5.88e+05, 6.10e+04,
+       7.71e+05, 1.87e+05, 8.14e+05, 2.81e+05, 1.34e+05, 4.20e+04,
+       6.80e+04, 5.20e+04, 4.50e+04, 4.10e+04])
+        self.f_res = np.array([6.290e+08, 8.400e+08, 1.066e+09, 1.076e+09, 1.608e+09, 1.884e+09,
+       2.218e+09, 2.533e+09, 2.547e+09, 2.782e+09, 3.008e+09, 3.223e+09,
+       3.284e+09, 3.463e+09, 3.643e+09, 3.761e+09, 3.900e+09, 4.000e+09,
+       4.080e+09, 4.210e+09, 1.076e+09, 1.100e+09, 1.955e+09, 2.075e+09,
+       2.118e+09, 2.199e+09, 2.576e+09, 2.751e+09, 3.370e+09, 5.817e+09,
+       5.817e+09, 1.210e+09, 1.280e+09, 1.415e+09, 1.415e+09, 1.415e+09,
+       1.415e+09, 1.395e+09, 1.401e+09, 1.570e+09, 1.610e+09, 1.620e+09,
+       1.861e+09, 1.890e+09, 2.495e+09, 6.960e+08, 9.100e+08, 1.069e+09,
+       1.078e+09, 1.155e+09, 1.232e+09, 1.343e+09])
+        self.Q_factor = np.array([ 500.,   10.,  500.,  500.,   40.,  500.,   15.,  384.,  340.,
+         20.,  450.,  512.,  600.,  805., 1040.,  965.,   50., 1300.,
+        600.,  200.,  700.,  700.,  450.,  600.,  600.,  750., 1000.,
+        500.,   30., 1000.,   10.,  315.,  200.,   75.,  270.,   75.,
+        270.,  200., 1100.,   55.,  980.,   60.,  810.,  175., 1190.,
+       7400., 8415., 7980., 7810., 6660., 5870., 7820.])
+        
+
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
 
@@ -244,7 +257,6 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
         import math 
-        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
         
         gamma_transition = 1 / np.sqrt(0.00192)  
         mom_compaction = 1 / gamma_transition ** 2
@@ -286,10 +298,30 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        table = np.loadtxt(this_directory + './EX_05_new_HQ_table.dat', comments='!')
-        self.R_shunt = table[:, 2] * 10 ** 6
-        self.f_res = table[:, 0] * 10 ** 9
-        self.Q_factor = table[:, 1]
+        self.R_shunt = np.array([3.88e+05, 6.00e+02, 6.98e+04, 8.72e+04, 2.04e+04, 1.09e+05,
+       9.10e+03, 1.80e+03, 7.50e+03, 2.60e+03, 2.80e+04, 2.90e+04,
+       1.20e+05, 7.90e+04, 4.80e+04, 1.45e+04, 5.00e+03, 1.05e+04,
+       1.25e+04, 1.00e+04, 4.30e+04, 4.50e+04, 1.15e+04, 4.50e+03,
+       6.90e+03, 1.54e+04, 1.10e+04, 4.70e+03, 8.50e+03, 1.70e+04,
+       5.00e+03, 6.33e+05, 4.99e+05, 3.64e+05, 1.77e+05, 6.14e+04,
+       2.99e+04, 3.79e+05, 2.43e+05, 1.74e+04, 5.88e+05, 6.10e+04,
+       7.71e+05, 1.87e+05, 8.14e+05, 2.81e+05, 1.34e+05, 4.20e+04,
+       6.80e+04, 5.20e+04, 4.50e+04, 4.10e+04])
+        self.f_res = np.array([6.290e+08, 8.400e+08, 1.066e+09, 1.076e+09, 1.608e+09, 1.884e+09,
+       2.218e+09, 2.533e+09, 2.547e+09, 2.782e+09, 3.008e+09, 3.223e+09,
+       3.284e+09, 3.463e+09, 3.643e+09, 3.761e+09, 3.900e+09, 4.000e+09,
+       4.080e+09, 4.210e+09, 1.076e+09, 1.100e+09, 1.955e+09, 2.075e+09,
+       2.118e+09, 2.199e+09, 2.576e+09, 2.751e+09, 3.370e+09, 5.817e+09,
+       5.817e+09, 1.210e+09, 1.280e+09, 1.415e+09, 1.415e+09, 1.415e+09,
+       1.415e+09, 1.395e+09, 1.401e+09, 1.570e+09, 1.610e+09, 1.620e+09,
+       1.861e+09, 1.890e+09, 2.495e+09, 6.960e+08, 9.100e+08, 1.069e+09,
+       1.078e+09, 1.155e+09, 1.232e+09, 1.343e+09])
+        self.Q_factor = np.array([ 500.,   10.,  500.,  500.,   40.,  500.,   15.,  384.,  340.,
+         20.,  450.,  512.,  600.,  805., 1040.,  965.,   50., 1300.,
+        600.,  200.,  700.,  700.,  450.,  600.,  600.,  750., 1000.,
+        500.,   30., 1000.,   10.,  315.,  200.,   75.,  270.,   75.,
+        270.,  200., 1100.,   55.,  980.,   60.,  810.,  175., 1190.,
+       7400., 8415., 7980., 7810., 6660., 5870., 7820.])
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
              
  
-- 
GitLab


From 74d0f41eb57d1f3227cc63485d261b365c255e74 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 15:44:16 +0100
Subject: [PATCH 11/48] test cases

---
 unittests/impedances/test_impedance.py | 8 +-------
 1 file changed, 1 insertion(+), 7 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 60b5c679..53b2c23a 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -122,21 +122,17 @@ class TestInducedVoltageTime(unittest.TestCase):
 class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
     def test_init(self):
-        # TODO Improve testcases
-    
+
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
 
     def test_init_mtw(self):
-        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       multi_turn_wake=True
                                       )
         
     def test_mtw_false_induced_volt(self):
-        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       multi_turn_wake=False)
-    
         ivr.process()
         ivr.induced_voltage_1turn()
         my_array = ivr.induced_voltage
@@ -144,7 +140,6 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
 
 
     def test_multi_rf_station(self):
-        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
                                       rf_station=self.rf_station, multi_turn_wake=False)
         ivr.process()
@@ -155,7 +150,6 @@ class  BaseTestInducedVoltageResonator(unittest.TestCase):
         
 
     def test_mtw_true_induced_voltage(self):
-        # TODO Improve testcases
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       rf_station=self.rf_station, multi_turn_wake=True)
         ivr.process()
-- 
GitLab


From 8a554a199b246126bb06f579ce9c813721cbbddc Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Mon, 16 Dec 2024 17:33:38 +0100
Subject: [PATCH 12/48] added mtw for resonator

---
 blond/impedances/impedance.py | 29 +++++++++++++++++++++++++++++
 1 file changed, 29 insertions(+)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 99deb89d..bfc74327 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -958,6 +958,7 @@ class InducedVoltageResonator(_InducedVoltage):
                  timeArray: Optional[NDArray] = None,
                  rf_station: Optional[RFStation] = None,
                  mtw_mode: Optional[MtwModeTypes] = None,
+                 array_length: Optional[int] = None,
                  use_regular_fft: bool = True) -> None:
 
         # Test if one or more quality factors is smaller than 0.5.
@@ -980,6 +981,7 @@ class InducedVoltageResonator(_InducedVoltage):
             self.tArray = timeArray
             self.atLineDensityTimes = False
 
+        self.array_length = array_length
         # Length of timeArray
         self.n_time = len(self.tArray)
 
@@ -1071,7 +1073,34 @@ class InducedVoltageResonator(_InducedVoltage):
                                 * self.beam.n_macroparticles * self.beam.ratio
         self.induced_voltage = (self.induced_voltage
                                 .astype(dtype=bm.precision.real_t, order='C', copy=False))
+        
 
+    def induced_voltage_mtw(self, beam_spectrum_dict: Optional[dict] = None):
+        """
+        Method to calculate the induced voltage taking into account the effect
+        from previous passages (multi-turn wake) @ F. Batsch
+        """
+        print('inside induced voltage mtw')
+        # work-around fix for timeshift, shift the entries in the array by 1 t_rev=512px
+        self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.profile.n_slices))
+        self.mtw_memory = self.mtw_memory[self.profile.n_slices:]
+
+        # Induced voltage of the current turn calculation
+        if beam_spectrum_dict is None:
+            beam_spectrum_dict = dict()
+            
+        self.induced_voltage_1turn(beam_spectrum_dict)
+
+        # Setting to zero to the last part to remove the contribution from the
+        # front wake
+        self.induced_voltage[self.n_induced_voltage - self.front_wake_buffer:] = 0
+
+        # Add the induced voltage of the current turn to the memory from previous
+        self.mtw_memory[:int(self.n_induced_voltage)] += self.induced_voltage  
+
+        self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
+        
+    
     # Implementation of Heaviside function
     def Heaviside(self, x):
         r"""
-- 
GitLab


From 492c4e767efbcc28bf1d4f828911a142eb82750b Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Wed, 18 Dec 2024 10:51:23 +0100
Subject: [PATCH 13/48] add mtw method for resonator

---
 blond/impedances/impedance.py | 39 ++++++++++++++++++-----------------
 1 file changed, 20 insertions(+), 19 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index bfc74327..da88b124 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -411,7 +411,10 @@ class _InducedVoltage:
         Method to calculate the induced voltage taking into account the effect
         from previous passages (multi-turn wake)
         """
+
+
         if beam_spectrum_dict is None:
+            print('beam_spectrum_dict is None')
             beam_spectrum_dict = dict()
         # Shift of the memory wake field by the current revolution period
         self.shift_trev()
@@ -419,17 +422,17 @@ class _InducedVoltage:
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
 
+
         # Setting to zero to the last part to remove the contribution from the
         # front wake
         self.induced_voltage[self.n_induced_voltage -
                              self.front_wake_buffer:] = 0
+        
 
         # Add the induced voltage of the current turn to the memory from
         # previous turns
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
-
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
-
     def shift_trev_freq(self):
         """
         Method to shift the induced voltage by a revolution period in the
@@ -452,12 +455,14 @@ class _InducedVoltage:
         """
 
         t_rev = self.rf_params.t_rev[self.rf_params.counter[0]]
-
+ 
         # self.mtw_memory = bm.interp_const_space(self.time_mtw + t_rev,
         self.mtw_memory = bm.interp(self.time_mtw + t_rev,
                                     self.time_mtw, self.mtw_memory,
                                     left=0, right=0)
 
+
+
     def _track(self):
         """
         Tracking method
@@ -1075,31 +1080,27 @@ class InducedVoltageResonator(_InducedVoltage):
                                 .astype(dtype=bm.precision.real_t, order='C', copy=False))
         
 
-    def induced_voltage_mtw(self, beam_spectrum_dict: Optional[dict] = None):
+    def induced_voltage_mtw(self, beam_spectrum_dict={}):
         """
         Method to calculate the induced voltage taking into account the effect
-        from previous passages (multi-turn wake) @ F. Batsch
+        from previous passages (multi-turn wake)
         """
-        print('inside induced voltage mtw')
-        # work-around fix for timeshift, shift the entries in the array by 1 t_rev=512px
-        self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.profile.n_slices))
-        self.mtw_memory = self.mtw_memory[self.profile.n_slices:]
+
+        # shift the entries in array by 1 t_rev=512px
+        self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.array_length))
+        self.mtw_memory = self.mtw_memory[self.array_length:]
 
         # Induced voltage of the current turn calculation
-        if beam_spectrum_dict is None:
-            beam_spectrum_dict = dict()
-            
         self.induced_voltage_1turn(beam_spectrum_dict)
 
-        # Setting to zero to the last part to remove the contribution from the
-        # front wake
-        self.induced_voltage[self.n_induced_voltage - self.front_wake_buffer:] = 0
+        # Add induced voltage of the current turn to the memory from previous
 
-        # Add the induced voltage of the current turn to the memory from previous
-        self.mtw_memory[:int(self.n_induced_voltage)] += self.induced_voltage  
-
-        self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
+        self.mtw_memory[:int(
+            self.n_time)] += self.induced_voltage  
+    
+        self.induced_voltage = self.mtw_memory[:self.n_time]
         
+
     
     # Implementation of Heaviside function
     def Heaviside(self, x):
-- 
GitLab


From 1bdd7d024f14b0ac23a39b52b46f136c5fa3df15 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Wed, 18 Dec 2024 14:30:54 +0100
Subject: [PATCH 14/48] remove unnecessary parameters from
 InducedVoltageResonator

---
 blond/impedances/impedance.py | 38 +++++++++++++++++------------------
 1 file changed, 19 insertions(+), 19 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index da88b124..99efbdca 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -911,7 +911,7 @@ class InducedVoltageResonator(_InducedVoltage):
     solving the convolution integral with the resonator impedance analytically.
     The line density need NOT be sampled at equidistant points. The times where
     the induced voltage is calculated need to be the same where the line
-    density is sampled. If no timeArray is passed, the induced voltage is
+    density is sampled. If no time_array is passed, the induced voltage is
     evaluated at the points of the line density. This is nececassry of
     compatability with other functions that calculate the induced voltage.
     Currently, it requires the all quality factors :math:`Q>0.5`
@@ -925,10 +925,13 @@ class InducedVoltageResonator(_InducedVoltage):
         Profile object
     resonators : Resonators
         Resonators object
-    timeArray : float array, optional
+    time_array : float array, optional
         Array of time values where the induced voltage is calculated.
         If left out, the induced voltage is calculated at the times of the line
         density.
+    array_length : int, optional
+        length of an array of one turn 
+    
 
     Attributes
     ----------
@@ -960,9 +963,8 @@ class InducedVoltageResonator(_InducedVoltage):
                  frequency_resolution: Optional[float] = None,
                  wake_length: Optional[float] = None,
                  multi_turn_wake: bool = False,
-                 timeArray: Optional[NDArray] = None,
+                 time_array: Optional[NDArray] = None,
                  rf_station: Optional[RFStation] = None,
-                 mtw_mode: Optional[MtwModeTypes] = None,
                  array_length: Optional[int] = None,
                  use_regular_fft: bool = True) -> None:
 
@@ -979,15 +981,15 @@ class InducedVoltageResonator(_InducedVoltage):
         # Optional array of time values where the induced voltage is calculated.
         # If left out, the induced voltage is calculated at the times of the
         # line density.
-        if timeArray is None:
+        if time_array is None:
             self.tArray = self.profile.bin_centers
             self.atLineDensityTimes = True
         else:
-            self.tArray = timeArray
+            self.tArray = time_array
             self.atLineDensityTimes = False
 
         self.array_length = array_length
-        # Length of timeArray
+
         self.n_time = len(self.tArray)
 
         # Copy of the shunt impedances of the Resonators in* :math:`\Omega`
@@ -1022,7 +1024,7 @@ class InducedVoltageResonator(_InducedVoltage):
                          frequency_resolution=frequency_resolution,
                          wake_length=wake_length,
                          multi_turn_wake=multi_turn_wake,
-                         rf_station=rf_station, mtw_mode=mtw_mode,
+                         rf_station=rf_station, mtw_mode=None,
                          use_regular_fft=use_regular_fft)
 
     def process(self):
@@ -1079,25 +1081,23 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = (self.induced_voltage
                                 .astype(dtype=bm.precision.real_t, order='C', copy=False))
         
-
+    
     def induced_voltage_mtw(self, beam_spectrum_dict={}):
+
         """
-        Method to calculate the induced voltage taking into account the effect
-        from previous passages (multi-turn wake)
+        Induced voltage method for InducedVoltageResonator. 
+        mtw_memory is shifted by one turn, setting the final values in the array as 0.
+        @fbatsch
         """
-
-        # shift the entries in array by 1 t_rev=512px
+        
+        # shift the entries in array by 1 t_rev and set to 0
         self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.array_length))
+        # remove one turn length of memory 
         self.mtw_memory = self.mtw_memory[self.array_length:]
-
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
-
         # Add induced voltage of the current turn to the memory from previous
-
-        self.mtw_memory[:int(
-            self.n_time)] += self.induced_voltage  
-    
+        self.mtw_memory[:int(self.n_time)] += self.induced_voltage  
         self.induced_voltage = self.mtw_memory[:self.n_time]
         
 
-- 
GitLab


From 00ab716a356191dbf646962d2c97bb4e203fda03 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Thu, 9 Jan 2025 16:57:59 +0100
Subject: [PATCH 15/48] edit tests

---
 unittests/impedances/test_impedance.py | 192 ++++++++++++-------------
 1 file changed, 93 insertions(+), 99 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 53b2c23a..9f53ea5d 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -16,7 +16,6 @@ Unittest for impedances.impedance
 import unittest
 
 import numpy as np
-
 from blond.beam.profile import CutOptions, Profile
 from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime, InducedVoltageResonator
 from blond.impedances.impedance_sources import Resonators
@@ -119,50 +118,10 @@ class TestInducedVoltageTime(unittest.TestCase):
         np.testing.assert_allclose(test_object.wake_length_input, 11e-9)
 
 
-class  BaseTestInducedVoltageResonator(unittest.TestCase):
-
-    def test_init(self):
-
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
-
-    def test_init_mtw(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=True
-                                      )
-        
-    def test_mtw_false_induced_volt(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=False)
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
-
-
-    def test_multi_rf_station(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
-                                      rf_station=self.rf_station, multi_turn_wake=False)
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
-
-        
-
-    def test_mtw_true_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station, multi_turn_wake=True)
-        ivr.process()
-        ivr.induced_voltage_mtw()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
-
-
-
-
-class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.TestCase):
+class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
 
     def setUp(self):
+        import os
         import numpy as np
         from blond.beam.beam import Beam, Proton
         from blond.beam.distributions import bigaussian
@@ -170,19 +129,18 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
         from blond.impedances.impedance_sources import Resonators
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
-        import math 
-        
+        import math
+
         gamma_transition = 1 / np.sqrt(0.00192)  
         mom_compaction = 1 / gamma_transition ** 2
         circ = 6911.56
-        
+
         ring = Ring(
             ring_length=circ,
             alpha_0=mom_compaction,
             synchronous_data=25.92e9,
             particle=Proton(),
             n_turns=2)
-                        
         
         self.rf_station = RFStation(
             ring=ring,
@@ -191,7 +149,7 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
             phi_rf_d=[0.0],
             n_rf=1,
         )
-        
+
         beam = Beam(
             ring=ring,
             n_macroparticles=1001,
@@ -204,42 +162,57 @@ class TestInducedVoltageResonatorRf1(BaseTestInducedVoltageResonator, unittest.T
             sigma_dt=2e-9 / 4, seed=1,
         )
         self.beam = beam
+
         cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        self.R_shunt = np.array([3.88e+05, 6.00e+02, 6.98e+04, 8.72e+04, 2.04e+04, 1.09e+05,
-       9.10e+03, 1.80e+03, 7.50e+03, 2.60e+03, 2.80e+04, 2.90e+04,
-       1.20e+05, 7.90e+04, 4.80e+04, 1.45e+04, 5.00e+03, 1.05e+04,
-       1.25e+04, 1.00e+04, 4.30e+04, 4.50e+04, 1.15e+04, 4.50e+03,
-       6.90e+03, 1.54e+04, 1.10e+04, 4.70e+03, 8.50e+03, 1.70e+04,
-       5.00e+03, 6.33e+05, 4.99e+05, 3.64e+05, 1.77e+05, 6.14e+04,
-       2.99e+04, 3.79e+05, 2.43e+05, 1.74e+04, 5.88e+05, 6.10e+04,
-       7.71e+05, 1.87e+05, 8.14e+05, 2.81e+05, 1.34e+05, 4.20e+04,
-       6.80e+04, 5.20e+04, 4.50e+04, 4.10e+04])
-        self.f_res = np.array([6.290e+08, 8.400e+08, 1.066e+09, 1.076e+09, 1.608e+09, 1.884e+09,
-       2.218e+09, 2.533e+09, 2.547e+09, 2.782e+09, 3.008e+09, 3.223e+09,
-       3.284e+09, 3.463e+09, 3.643e+09, 3.761e+09, 3.900e+09, 4.000e+09,
-       4.080e+09, 4.210e+09, 1.076e+09, 1.100e+09, 1.955e+09, 2.075e+09,
-       2.118e+09, 2.199e+09, 2.576e+09, 2.751e+09, 3.370e+09, 5.817e+09,
-       5.817e+09, 1.210e+09, 1.280e+09, 1.415e+09, 1.415e+09, 1.415e+09,
-       1.415e+09, 1.395e+09, 1.401e+09, 1.570e+09, 1.610e+09, 1.620e+09,
-       1.861e+09, 1.890e+09, 2.495e+09, 6.960e+08, 9.100e+08, 1.069e+09,
-       1.078e+09, 1.155e+09, 1.232e+09, 1.343e+09])
-        self.Q_factor = np.array([ 500.,   10.,  500.,  500.,   40.,  500.,   15.,  384.,  340.,
-         20.,  450.,  512.,  600.,  805., 1040.,  965.,   50., 1300.,
-        600.,  200.,  700.,  700.,  450.,  600.,  600.,  750., 1000.,
-        500.,   30., 1000.,   10.,  315.,  200.,   75.,  270.,   75.,
-        270.,  200., 1100.,   55.,  980.,   60.,  810.,  175., 1190.,
-       7400., 8415., 7980., 7810., 6660., 5870., 7820.])
         
+        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+        table = np.loadtxt(this_directory + './data/example_data.dat', comments='!')
+        self.R_shunt = table[:, 2] * 10**6
+        self.f_res = table[:, 0] * 10**9
+        self.Q_factor = table[:, 1]
 
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
 
+    def test_init(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
 
+    def test_init_mtw(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=True
+                                      )
+        
+    def test_mtw_false_induced_volt(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=False)
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+    def test_multi_rf_station(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
+                                      rf_station=self.rf_station, multi_turn_wake=False)
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
 
-class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.TestCase):
+    def test_mtw_true_induced_voltage(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      rf_station=self.rf_station, multi_turn_wake=True)
+        ivr.process()
+        ivr.induced_voltage_mtw()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+
+
+
+class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
 
     def setUp(self):
         import os
@@ -292,33 +265,54 @@ class TestInducedVoltageResonatorRf2(BaseTestInducedVoltageResonator, unittest.T
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        self.R_shunt = np.array([3.88e+05, 6.00e+02, 6.98e+04, 8.72e+04, 2.04e+04, 1.09e+05,
-       9.10e+03, 1.80e+03, 7.50e+03, 2.60e+03, 2.80e+04, 2.90e+04,
-       1.20e+05, 7.90e+04, 4.80e+04, 1.45e+04, 5.00e+03, 1.05e+04,
-       1.25e+04, 1.00e+04, 4.30e+04, 4.50e+04, 1.15e+04, 4.50e+03,
-       6.90e+03, 1.54e+04, 1.10e+04, 4.70e+03, 8.50e+03, 1.70e+04,
-       5.00e+03, 6.33e+05, 4.99e+05, 3.64e+05, 1.77e+05, 6.14e+04,
-       2.99e+04, 3.79e+05, 2.43e+05, 1.74e+04, 5.88e+05, 6.10e+04,
-       7.71e+05, 1.87e+05, 8.14e+05, 2.81e+05, 1.34e+05, 4.20e+04,
-       6.80e+04, 5.20e+04, 4.50e+04, 4.10e+04])
-        self.f_res = np.array([6.290e+08, 8.400e+08, 1.066e+09, 1.076e+09, 1.608e+09, 1.884e+09,
-       2.218e+09, 2.533e+09, 2.547e+09, 2.782e+09, 3.008e+09, 3.223e+09,
-       3.284e+09, 3.463e+09, 3.643e+09, 3.761e+09, 3.900e+09, 4.000e+09,
-       4.080e+09, 4.210e+09, 1.076e+09, 1.100e+09, 1.955e+09, 2.075e+09,
-       2.118e+09, 2.199e+09, 2.576e+09, 2.751e+09, 3.370e+09, 5.817e+09,
-       5.817e+09, 1.210e+09, 1.280e+09, 1.415e+09, 1.415e+09, 1.415e+09,
-       1.415e+09, 1.395e+09, 1.401e+09, 1.570e+09, 1.610e+09, 1.620e+09,
-       1.861e+09, 1.890e+09, 2.495e+09, 6.960e+08, 9.100e+08, 1.069e+09,
-       1.078e+09, 1.155e+09, 1.232e+09, 1.343e+09])
-        self.Q_factor = np.array([ 500.,   10.,  500.,  500.,   40.,  500.,   15.,  384.,  340.,
-         20.,  450.,  512.,  600.,  805., 1040.,  965.,   50., 1300.,
-        600.,  200.,  700.,  700.,  450.,  600.,  600.,  750., 1000.,
-        500.,   30., 1000.,   10.,  315.,  200.,   75.,  270.,   75.,
-        270.,  200., 1100.,   55.,  980.,   60.,  810.,  175., 1190.,
-       7400., 8415., 7980., 7810., 6660., 5870., 7820.])
+        
+        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
+        table = np.loadtxt(this_directory + './data/example_data.dat', comments='!')
+        self.R_shunt = table[:, 2] * 10**6
+        self.f_res = table[:, 0] * 10
+        self.Q_factor = table[:, 1]
+
+
+
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
-             
- 
+
+
+    def test_init(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
+
+    def test_init_mtw(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=True
+                                      )
+        
+    def test_mtw_false_induced_volt(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      multi_turn_wake=False)
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+    def test_multi_rf_station(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
+                                      rf_station=self.rf_station, multi_turn_wake=False)
+        ivr.process()
+        ivr.induced_voltage_1turn()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+    def test_mtw_true_induced_voltage(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      rf_station=self.rf_station, multi_turn_wake=True)
+        ivr.process()
+        ivr.induced_voltage_mtw()
+        my_array = ivr.induced_voltage
+        self.assertTrue(np.any(my_array != 0.0))
+
+
+
+
+
 
 
 if __name__ == '__main__':
-- 
GitLab


From 715d9156bcbbabaeeb69021366a2e376fa48d693 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Thu, 9 Jan 2025 17:08:52 +0100
Subject: [PATCH 16/48] add test example data

---
 unittests/impedances/data/Q_factor.npy | Bin 0 -> 544 bytes
 unittests/impedances/data/R_shunt.npy  | Bin 0 -> 544 bytes
 unittests/impedances/data/f_res.npy    | Bin 0 -> 544 bytes
 3 files changed, 0 insertions(+), 0 deletions(-)
 create mode 100644 unittests/impedances/data/Q_factor.npy
 create mode 100644 unittests/impedances/data/R_shunt.npy
 create mode 100644 unittests/impedances/data/f_res.npy

diff --git a/unittests/impedances/data/Q_factor.npy b/unittests/impedances/data/Q_factor.npy
new file mode 100644
index 0000000000000000000000000000000000000000..1132114819c01ac1a81b9235d42d363a97760e18
GIT binary patch
literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<&-aHw|x(+nyQ8pemw3@)hR40=%Y6;Sa~DBlD^E7U;w
z4G{W36NJ|2gwPHXp!_}v&EN^41ExS|n0kgx2>(GVDh)FqrVd6Iv_sT6^h0TwzwA)O
zVgB4u3=xOfb07%9e^3CS(fKfS1`{Fj4c-uX0n8tM5PksEpAS+Y{27xV^ou=Ux*`7n
Ygs-z7LO1P$&@8(k^t){kx@VsQ02D)1r2qf`

literal 0
HcmV?d00001

diff --git a/unittests/impedances/data/R_shunt.npy b/unittests/impedances/data/R_shunt.npy
new file mode 100644
index 0000000000000000000000000000000000000000..f45cbe077aab4a0896cd1cd41bde6b2ca7ccbd04
GIT binary patch
literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<(2Xjmui2%;GdG&z8%2Ck1_n!)QUgnn}wOgGH=1>qk%
z2%#0`Kxm)65ISH9gpRlkp@Z%~<^O`|2G7qBe!+97{22({u^B=)9f8nUC!peoA@u(z
z5c=m+2wipzLQmNUp}*{c&?4s`bj?u+9kU5S&pH65FF<IheF1^WV1GX-kpt5W=3-#l
zVKyt6X88RUqTcu}n0}xL38w=qc_H%uAm$zTs0iVogs3lQR|oSCTxJ8)3D-5i^a3?O
cFx?Q%0Hzt5p#GZm0b)<eO9&0~C+iaj03BFuSpWb4

literal 0
HcmV?d00001

diff --git a/unittests/impedances/data/f_res.npy b/unittests/impedances/data/f_res.npy
new file mode 100644
index 0000000000000000000000000000000000000000..72fff53c63ff2c43e265be3a912da52b2d4479b2
GIT binary patch
literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<&zh6RuJ9dZQH3>!I4f@p{LN6v%kxf~ZjbOXn?>md5T
zW1%}B`a*r%0}!1NaQ_jAj$ru(rVH$LJpu6<{4$?`=!zuf7a)4U?s72wL-p)S5PwF6
z%4-mPqLcLvh<>nw?=6UCNIeUt8-Az01Mvku#J&g74p4U-&`JTjXTxsAyC89c`q%eB
zbV7UoeGpyYlmvEHLx|yH5PyekJ=i@L4oidWWocXe6(kO$7o1?f1mZJX5V{PaKb#i5
zijUsVaPkUBULfTGm{t%@1&7Omm4Css!+GT!Ao&YwyKjSNh2xsw@Y->88aR9y`1c(F
mi67W6cN#=5xO5vF&Ixt`5P$tNy$BL#*cuIXSAjC~6-NMe6rXkg

literal 0
HcmV?d00001

-- 
GitLab


From d1baa3535cd4b3f00481f4a490780cfcdebf8069 Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Thu, 9 Jan 2025 17:10:31 +0100
Subject: [PATCH 17/48] read test data in test impedance

---
 unittests/impedances/test_impedance.py | 7 +++----
 1 file changed, 3 insertions(+), 4 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 9f53ea5d..ab2a0f73 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -169,10 +169,9 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
                                FitOptions(fit_option='gaussian'))
         
         this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        table = np.loadtxt(this_directory + './data/example_data.dat', comments='!')
-        self.R_shunt = table[:, 2] * 10**6
-        self.f_res = table[:, 0] * 10**9
-        self.Q_factor = table[:, 1]
+        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10**6
+        self.f_res = np.load(this_directory + './data/f_res.npy') * 10**9
+        self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
 
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
-- 
GitLab


From 7d8f9695d495eac951a6d44a2c326ec8cfd06e7a Mon Sep 17 00:00:00 2001
From: Elleanor Rose Lamb <elamb@pcsy107.dyndns.cern.ch>
Date: Fri, 10 Jan 2025 15:09:52 +0100
Subject: [PATCH 18/48] edit tests

---
 unittests/impedances/test_impedance.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index ab2a0f73..f081af80 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -266,10 +266,10 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
                                FitOptions(fit_option='gaussian'))
         
         this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        table = np.loadtxt(this_directory + './data/example_data.dat', comments='!')
-        self.R_shunt = table[:, 2] * 10**6
-        self.f_res = table[:, 0] * 10
-        self.Q_factor = table[:, 1]
+        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10**6
+        self.f_res = np.load(this_directory + './data/f_res.npy') * 10**9
+        self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
+
 
 
 
-- 
GitLab


From 6bbaa66f5e0c6d19852154c19a95079cea17c801 Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Fri, 10 Jan 2025 17:16:09 +0100
Subject: [PATCH 19/48] Added legacy_support for timeArray

---
 blond/utils/legacy_support.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/blond/utils/legacy_support.py b/blond/utils/legacy_support.py
index b92bd5ef..bab98d4b 100644
--- a/blond/utils/legacy_support.py
+++ b/blond/utils/legacy_support.py
@@ -73,6 +73,7 @@ __new_by_old = {
     "RingAndRFSection_list": "ring_and_rf_section",
     "BeamFeedback": "beam_feedback",
     "CavityFeedback": "cavity_feedback",
-    "smoothOption": "smooth_option"
+    "smoothOption": "smooth_option",
+    "timeArray": "time_array",
 }
 handle_legacy_kwargs = _handle_legacy_kwargs(__new_by_old)
-- 
GitLab


From b5bdd315faa5bcc568ad2ed831e8b66e94c5d4ad Mon Sep 17 00:00:00 2001
From: slauber <simon.fabian.lauber@cern.ch>
Date: Fri, 10 Jan 2025 10:26:26 +0100
Subject: [PATCH 20/48] Deprecated C++ cumtrapz implementation, as it doesnt
 match scipy.integrate.cumulative_trapezoid anymore. Altered comparison of
 expected/actual values in three tests to make them working again.

---
 blond/compile.py                              |  3 ++-
 blond/utils/bmath.py                          |  2 +-
 blond/utils/butils_wrap_cpp.py                |  3 +++
 .../beam_profile/test_beam_profile_object.py  |  4 +---
 unittests/llrf/test_cavity_feedback.py        | 24 +++++++++----------
 unittests/utils/test_blondmath.py             |  8 +++----
 6 files changed, 23 insertions(+), 21 deletions(-)

diff --git a/blond/compile.py b/blond/compile.py
index 5314d76d..069dd1ea 100644
--- a/blond/compile.py
+++ b/blond/compile.py
@@ -310,9 +310,10 @@ def compile_cuda_library(args, nvccflags, float_flags, cuda_files, nvcc):
     # Compile the GPU library
     # print('\n' + ''.join(['='] * 80))
     print('\nCompiling the CUDA library')
+    import cupy as cp
+
     if args['gpu'] == 'discover':
         print('Discovering the device compute capability..')
-        import cupy as cp
 
         dev = cp.cuda.Device(0)
         dev_name = cp.cuda.runtime.getDeviceProperties(dev)['name']
diff --git a/blond/utils/bmath.py b/blond/utils/bmath.py
index 636dbf6c..f0a44623 100644
--- a/blond/utils/bmath.py
+++ b/blond/utils/bmath.py
@@ -47,7 +47,7 @@ def use_cpp() -> None:
         'interp_cpp': _cpp.interp_cpp,
         'interp_const_space': _cpp.interp_const_space,
         'interp_const_bin': _cpp.interp_const_bin,
-        'cumtrapz': _cpp.cumtrapz,
+        # 'cumtrapz': _cpp.cumtrapz, # needs reimplementation to match cupy/scipy
         'trapz_cpp': _cpp.trapz_cpp,
         'linspace_cpp': _cpp.linspace_cpp,
         'argmin_cpp': _cpp.argmin_cpp,
diff --git a/blond/utils/butils_wrap_cpp.py b/blond/utils/butils_wrap_cpp.py
index 6e658789..9e5ea784 100644
--- a/blond/utils/butils_wrap_cpp.py
+++ b/blond/utils/butils_wrap_cpp.py
@@ -530,6 +530,9 @@ def irfft_packed(signal: NDArray, fftsize: int = 0, result: Optional[NDArray] =
 
 def cumtrapz(y: NDArray, x: Optional[NDArray] = None, dx: float = 1.0, initial: Optional[float] = None,
              result: Optional[NDArray] = None) -> NDArray:
+    raise NotImplementedError("The bmath.cumtrapz behaviour differs from"
+                              "scipy.integrate.cumulative_trapezoid."
+                              "Please contact the developers if you need this routine!")
     if x is not None:
         # IntegrationError
         raise RuntimeError('[cumtrapz] x attribute is not yet supported')
diff --git a/unittests/beam_profile/test_beam_profile_object.py b/unittests/beam_profile/test_beam_profile_object.py
index e9cd4925..8c084b0b 100644
--- a/unittests/beam_profile/test_beam_profile_object.py
+++ b/unittests/beam_profile/test_beam_profile_object.py
@@ -208,7 +208,7 @@ class testProfileClass(unittest.TestCase):
              0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
              0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
              0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.823322140525,
-             0.176677859475, 0.0, -0.255750723657, 1.58255526749, 5.591836814,
+             0.17667805113, 0.0, -0.255750723657, 1.58255526749, 5.591836814,
              8.52273531218, 26.7796551186, 34.3314288768, 51.0732749092,
              96.1902354197, 157.367580381, 225.376283353, 374.772960616,
              486.894110696, 747.13998452, 949.483971664, 1340.35510048,
@@ -278,8 +278,6 @@ class testProfileClass(unittest.TestCase):
             [9.27853156526e-08, 9.24434506817e-08, 9.18544356769e-08],
             rtol=rtol, atol=atol,
             err_msg='Bunch length values not correct')
-
-
 if __name__ == '__main__':
 
     unittest.main()
diff --git a/unittests/llrf/test_cavity_feedback.py b/unittests/llrf/test_cavity_feedback.py
index 82f2d322..c0d53593 100644
--- a/unittests/llrf/test_cavity_feedback.py
+++ b/unittests/llrf/test_cavity_feedback.py
@@ -444,16 +444,16 @@ class TestSPSOneTurnFeedback(unittest.TestCase):
         ref_DV_MOD_FR = np.load(os.path.join(this_directory, "ref_DV_MOD_FR.npy"))
 
         # Test real part
-        np.testing.assert_allclose(np.around(self.OTFB.DV_MOD_FR[-self.OTFB.n_coarse:].real, 12),
-                                   np.around(ref_DV_MOD_FR.real, 12),
-                                   rtol=1e-6, atol=0,
+        np.testing.assert_allclose(self.OTFB.DV_MOD_FR[-self.OTFB.n_coarse:].real,
+                                   ref_DV_MOD_FR.real,
+                                   rtol=1e-6, atol=1e-9,
                                    err_msg="In TestSPSOneTurnFeedback test_mod_to_fr(), "
                                            "mismatch in real part of modulated signal")
 
         # Test imaginary part
-        np.testing.assert_allclose(np.around(self.OTFB.DV_MOD_FR[-self.OTFB.n_coarse:].imag, 12),
-                                   np.around(ref_DV_MOD_FR.imag, 12),
-                                   rtol=1e-6, atol=0,
+        np.testing.assert_allclose(self.OTFB.DV_MOD_FR[-self.OTFB.n_coarse:].imag,
+                                   ref_DV_MOD_FR.imag,
+                                   rtol=1e-6, atol=1e-9,
                                    err_msg="In TestSPSOneTurnFeedback test_mod_to_fr(), "
                                            "mismatch in imaginary part of modulated signal")
 
@@ -495,16 +495,16 @@ class TestSPSOneTurnFeedback(unittest.TestCase):
         ref_DV_MOD_FRF = np.load(os.path.join(this_directory, "ref_DV_MOD_FRF.npy"))
 
         # Test real part
-        np.testing.assert_allclose(np.around(self.OTFB.DV_MOD_FRF[-self.OTFB.n_coarse:].real, 12),
-                                   np.around(ref_DV_MOD_FRF.real, 12),
-                                   rtol=1e-6, atol=0,
+        np.testing.assert_allclose(self.OTFB.DV_MOD_FRF[-self.OTFB.n_coarse:].real,
+                                   ref_DV_MOD_FRF.real,
+                                   rtol=1e-6, atol=1e-9,
                                    err_msg="In TestSPSOneTurnFeedback test_mod_to_frf(), "
                                            "mismatch in real part of modulated signal")
 
         # Test imaginary part
-        np.testing.assert_allclose(np.around(self.OTFB.DV_MOD_FRF[-self.OTFB.n_coarse:].imag, 12),
-                                   np.around(ref_DV_MOD_FRF.imag, 12),
-                                   rtol=1e-6, atol=0,
+        np.testing.assert_allclose(self.OTFB.DV_MOD_FRF[-self.OTFB.n_coarse:].imag,
+                                   ref_DV_MOD_FRF.imag,
+                                   rtol=1e-6, atol=1e-9,
                                    err_msg="In TestSPSOneTurnFeedback test_mod_to_frf(), "
                                            "mismatch in imaginary part of modulated signal")
 
diff --git a/unittests/utils/test_blondmath.py b/unittests/utils/test_blondmath.py
index a182fa63..7c777dc2 100644
--- a/unittests/utils/test_blondmath.py
+++ b/unittests/utils/test_blondmath.py
@@ -637,7 +637,7 @@ class TestCumTrapz(unittest.TestCase):
         pass
 
     def test_cumtrapz_1(self):
-        import scipy.integrate
+        self.skipTest("scipy.integrate.cumulative_trapezoid behaviour differs from cumtrapz")
         y = np.random.randn(100)
         initial = np.random.rand()
         np.testing.assert_almost_equal(bm.cumtrapz(y, initial=initial),
@@ -646,14 +646,14 @@ class TestCumTrapz(unittest.TestCase):
                                        decimal=8)
 
     def test_cumtrapz_2(self):
-        import scipy.integrate
+        self.skipTest("scipy.integrate.cumulative_trapezoid behaviour differs from cumtrapz")
         y = np.random.randn(100)
         np.testing.assert_almost_equal(bm.cumtrapz(y),
                                        cumtrapz(y),
                                        decimal=8)
 
     def test_cumtrapz_3(self):
-        import scipy.integrate
+        self.skipTest("scipy.integrate.cumulative_trapezoid behaviour differs from cumtrapz")
         y = np.random.randn(100)
         dx = np.random.rand()
         np.testing.assert_almost_equal(bm.cumtrapz(y, dx=dx),
@@ -661,7 +661,7 @@ class TestCumTrapz(unittest.TestCase):
                                        decimal=8)
 
     def test_cumtrapz_4(self):
-        import scipy.integrate
+        self.skipTest("scipy.integrate.cumulative_trapezoid behaviour differs from cumtrapz")
         y = np.random.randn(100)
         dx = np.random.rand()
         initial = np.random.rand()
-- 
GitLab


From fe3371f7b4bac9df604beb468be2176504540d1e Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Mon, 13 Jan 2025 13:47:50 +0100
Subject: [PATCH 21/48] Added test for multi station ring energy calculation

---
 blond/input_parameters/ring.py          |  4 +++-
 unittests/input_parameters/test_ring.py | 28 +++++++++++++++++++++++++
 2 files changed, 31 insertions(+), 1 deletion(-)

diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 57fb16a5..7aae56c9 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -282,10 +282,11 @@ class Ring:
         self.f_rev: NDArray = 1 / self.t_rev
         self.omega_rev: NDArray = 2 * np.pi * self.f_rev
 
+        
 
-                
         if self.n_sections == 1:
             self.delta_E = np.diff(self.energy, axis=1)
+
         else:
             self.delta_E = np.reshape(np.zeros(n_turns*n_sections),(n_sections,n_turns))
             for i in range(n_turns):  
@@ -296,6 +297,7 @@ class Ring:
                         self.delta_E[j,i] = self.energy[j,i+1]-self.energy[j-1,i+1]
 
 
+
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
             interp_time = self.cycle_time
diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index 71863da9..a9c33888 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -192,5 +192,33 @@ class TestGeneralParameters(unittest.TestCase):
         Ring(C, alpha, momentum, Proton(), n_turns)  # might crash
 
 
+
+    def test_bug_multi_RF_dE(self):
+        from blond.beam.beam import Proton
+        C = 2 * np.pi * 1100.009  # Ring circumference [m]
+        gamma_t = 18.0  # Gamma at transition
+        alpha = 1 / gamma_t ** 2  # Momentum compaction factor
+        n_turns = self.n_turns
+        n_sections = 8
+        l_per_section = C/ n_sections
+        section_lengths = np.full(n_sections, l_per_section)
+        mon_inj = 1*1e10
+        mon_ext = 10*1e10
+        momentum = np.zeros(n_sections * (n_turns + 1))
+        momentum[0:(n_sections - 1)] = mon_inj
+        # linear momentum program
+        momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, 
+                                                num=(n_sections * (n_turns + 1)) - (n_sections - 1))
+        momentum = momentum.reshape(n_turns + 1, n_sections).T
+        ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
+            synchronous_data=momentum, particle=Proton(), 
+            n_turns=n_turns, n_sections=n_sections,
+                    synchronous_data_type='momentum')
+        assert np.isclose(ring.delta_E[0, 0] / ring.delta_E[0, -1], 1,0.1), \
+        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[0, 0]} / {ring.delta_E[0, -1]}"
+
+
+
+
 if __name__ == '__main__':
     unittest.main()
-- 
GitLab


From 126df768fe18de46ef6146f2f3c97aa421cccccd Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Mon, 13 Jan 2025 16:41:23 +0100
Subject: [PATCH 22/48] added test for ring

---
 blond/input_parameters/ring.py          |  4 ++-
 elle_scripts/testing_ring.py            | 42 +++++++++++++++++++++++++
 unittests/input_parameters/test_ring.py | 21 ++++++-------
 3 files changed, 54 insertions(+), 13 deletions(-)
 create mode 100644 elle_scripts/testing_ring.py

diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 7aae56c9..96da0518 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -282,13 +282,13 @@ class Ring:
         self.f_rev: NDArray = 1 / self.t_rev
         self.omega_rev: NDArray = 2 * np.pi * self.f_rev
 
-        
 
         if self.n_sections == 1:
             self.delta_E = np.diff(self.energy, axis=1)
 
         else:
             self.delta_E = np.reshape(np.zeros(n_turns*n_sections),(n_sections,n_turns))
+            # to take the section from the previous turn
             for i in range(n_turns):  
                 for j in range(n_sections):
                     if j == 0:
@@ -298,6 +298,8 @@ class Ring:
 
 
 
+
+
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
             interp_time = self.cycle_time
diff --git a/elle_scripts/testing_ring.py b/elle_scripts/testing_ring.py
new file mode 100644
index 00000000..358586e6
--- /dev/null
+++ b/elle_scripts/testing_ring.py
@@ -0,0 +1,42 @@
+# %%
+
+from blond.beam.beam import Proton
+from matplotlib import pyplot as plt
+
+import sys
+import unittest
+
+import numpy as np
+
+from blond.beam.beam import Electron
+from blond.input_parameters.ring import Ring
+from blond.input_parameters.ring_options import convert_data
+
+C = 2 * np.pi * 1100.009  # Ring circumference [m]
+gamma_t = 18.0  # Gamma at transition
+alpha = 1 / gamma_t ** 2  # Momentum compaction factor
+n_turns = 10
+n_sections = 8
+l_per_section = C/ n_sections
+section_lengths = np.full(n_sections, l_per_section)
+mon_inj = 1*1e10
+mon_ext = 10*1e10
+
+momentum = np.zeros(n_sections * (n_turns + 1))
+momentum[0:(n_sections - 1)] = mon_inj
+# linear momentum program
+momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, num=(n_sections * (n_turns + 1)) - (n_sections - 1))
+momentum = momentum.reshape(n_turns + 1, n_sections).T
+
+
+# %%
+
+ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
+            synchronous_data=momentum, particle=Proton(), 
+            n_turns=n_turns, n_sections=n_sections,
+                    synchronous_data_type='momentum')
+
+
+print(ring.delta_E)
+
+# %%
diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index a9c33888..a2a12880 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -199,23 +199,20 @@ class TestGeneralParameters(unittest.TestCase):
         gamma_t = 18.0  # Gamma at transition
         alpha = 1 / gamma_t ** 2  # Momentum compaction factor
         n_turns = self.n_turns
-        n_sections = 8
+        n_sections = 4
         l_per_section = C/ n_sections
         section_lengths = np.full(n_sections, l_per_section)
-        mon_inj = 1*1e10
-        mon_ext = 10*1e10
-        momentum = np.zeros(n_sections * (n_turns + 1))
-        momentum[0:(n_sections - 1)] = mon_inj
-        # linear momentum program
-        momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, 
-                                                num=(n_sections * (n_turns + 1)) - (n_sections - 1))
-        momentum = momentum.reshape(n_turns + 1, n_sections).T
-        ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
+        # should be a linear energy ramp, taking the value from the section in the previous turn
+        momentum = np.arange(n_sections * (n_turns + 1),dtype=float).reshape(n_turns + 1, n_sections).T*1e10
+        momentum[:, 0] = momentum[0, 1]
+        ring = Ring(ring_length=section_lengths, alpha_0=alpha,
             synchronous_data=momentum, particle=Proton(), 
             n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
-        assert np.isclose(ring.delta_E[0, 0] / ring.delta_E[0, -1], 1,0.1), \
-        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[0, 0]} / {ring.delta_E[0, -1]}"
+        a = ring.delta_E[1,0]
+        b = ring.delta_E[1,1]
+        assert np.isclose(a/b,1,rtol=1e-1), \
+        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[1, 0]} == {ring.delta_E[1,1]}"
 
 
 
-- 
GitLab


From 1a026ced2b4c3f1ac73627d8e75af77f311b7979 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Mon, 13 Jan 2025 16:43:04 +0100
Subject: [PATCH 23/48] Revert "added test for ring"

This reverts commit 126df768fe18de46ef6146f2f3c97aa421cccccd.
---
 blond/input_parameters/ring.py          |  4 +--
 elle_scripts/testing_ring.py            | 42 -------------------------
 unittests/input_parameters/test_ring.py | 21 +++++++------
 3 files changed, 13 insertions(+), 54 deletions(-)
 delete mode 100644 elle_scripts/testing_ring.py

diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 96da0518..7aae56c9 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -282,13 +282,13 @@ class Ring:
         self.f_rev: NDArray = 1 / self.t_rev
         self.omega_rev: NDArray = 2 * np.pi * self.f_rev
 
+        
 
         if self.n_sections == 1:
             self.delta_E = np.diff(self.energy, axis=1)
 
         else:
             self.delta_E = np.reshape(np.zeros(n_turns*n_sections),(n_sections,n_turns))
-            # to take the section from the previous turn
             for i in range(n_turns):  
                 for j in range(n_sections):
                     if j == 0:
@@ -298,8 +298,6 @@ class Ring:
 
 
 
-
-
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
             interp_time = self.cycle_time
diff --git a/elle_scripts/testing_ring.py b/elle_scripts/testing_ring.py
deleted file mode 100644
index 358586e6..00000000
--- a/elle_scripts/testing_ring.py
+++ /dev/null
@@ -1,42 +0,0 @@
-# %%
-
-from blond.beam.beam import Proton
-from matplotlib import pyplot as plt
-
-import sys
-import unittest
-
-import numpy as np
-
-from blond.beam.beam import Electron
-from blond.input_parameters.ring import Ring
-from blond.input_parameters.ring_options import convert_data
-
-C = 2 * np.pi * 1100.009  # Ring circumference [m]
-gamma_t = 18.0  # Gamma at transition
-alpha = 1 / gamma_t ** 2  # Momentum compaction factor
-n_turns = 10
-n_sections = 8
-l_per_section = C/ n_sections
-section_lengths = np.full(n_sections, l_per_section)
-mon_inj = 1*1e10
-mon_ext = 10*1e10
-
-momentum = np.zeros(n_sections * (n_turns + 1))
-momentum[0:(n_sections - 1)] = mon_inj
-# linear momentum program
-momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, num=(n_sections * (n_turns + 1)) - (n_sections - 1))
-momentum = momentum.reshape(n_turns + 1, n_sections).T
-
-
-# %%
-
-ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
-            synchronous_data=momentum, particle=Proton(), 
-            n_turns=n_turns, n_sections=n_sections,
-                    synchronous_data_type='momentum')
-
-
-print(ring.delta_E)
-
-# %%
diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index a2a12880..a9c33888 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -199,20 +199,23 @@ class TestGeneralParameters(unittest.TestCase):
         gamma_t = 18.0  # Gamma at transition
         alpha = 1 / gamma_t ** 2  # Momentum compaction factor
         n_turns = self.n_turns
-        n_sections = 4
+        n_sections = 8
         l_per_section = C/ n_sections
         section_lengths = np.full(n_sections, l_per_section)
-        # should be a linear energy ramp, taking the value from the section in the previous turn
-        momentum = np.arange(n_sections * (n_turns + 1),dtype=float).reshape(n_turns + 1, n_sections).T*1e10
-        momentum[:, 0] = momentum[0, 1]
-        ring = Ring(ring_length=section_lengths, alpha_0=alpha,
+        mon_inj = 1*1e10
+        mon_ext = 10*1e10
+        momentum = np.zeros(n_sections * (n_turns + 1))
+        momentum[0:(n_sections - 1)] = mon_inj
+        # linear momentum program
+        momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, 
+                                                num=(n_sections * (n_turns + 1)) - (n_sections - 1))
+        momentum = momentum.reshape(n_turns + 1, n_sections).T
+        ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
             synchronous_data=momentum, particle=Proton(), 
             n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
-        a = ring.delta_E[1,0]
-        b = ring.delta_E[1,1]
-        assert np.isclose(a/b,1,rtol=1e-1), \
-        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[1, 0]} == {ring.delta_E[1,1]}"
+        assert np.isclose(ring.delta_E[0, 0] / ring.delta_E[0, -1], 1,0.1), \
+        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[0, 0]} / {ring.delta_E[0, -1]}"
 
 
 
-- 
GitLab


From 75b678b8fc90656227ca3d1807fa7cbf2bd87dc3 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Mon, 13 Jan 2025 16:46:42 +0100
Subject: [PATCH 24/48] added test for ring

---
 unittests/input_parameters/test_ring.py | 25 ++++++++++++-------------
 1 file changed, 12 insertions(+), 13 deletions(-)

diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index a9c33888..197d79ae 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -193,29 +193,28 @@ class TestGeneralParameters(unittest.TestCase):
 
 
 
+
     def test_bug_multi_RF_dE(self):
         from blond.beam.beam import Proton
         C = 2 * np.pi * 1100.009  # Ring circumference [m]
         gamma_t = 18.0  # Gamma at transition
         alpha = 1 / gamma_t ** 2  # Momentum compaction factor
         n_turns = self.n_turns
-        n_sections = 8
+        n_sections = 4
         l_per_section = C/ n_sections
         section_lengths = np.full(n_sections, l_per_section)
-        mon_inj = 1*1e10
-        mon_ext = 10*1e10
-        momentum = np.zeros(n_sections * (n_turns + 1))
-        momentum[0:(n_sections - 1)] = mon_inj
-        # linear momentum program
-        momentum[n_sections - 1:] = np.linspace(mon_inj, mon_ext, 
-                                                num=(n_sections * (n_turns + 1)) - (n_sections - 1))
-        momentum = momentum.reshape(n_turns + 1, n_sections).T
-        ring = Ring(ring_length=section_lengths, alpha_0=alpha, 
-            synchronous_data=momentum, particle=Proton(), 
+        # should be a linear energy ramp, taking the value from the section in the previous turn
+        momentum = np.arange(n_sections * (n_turns + 1),dtype=float).reshape(n_turns + 1, n_sections).T*1e10
+        momentum[:, 0] = momentum[0, 1]
+        ring = Ring(ring_length=section_lengths, alpha_0=alpha,
+            synchronous_data=momentum, particle=Proton(),
             n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
-        assert np.isclose(ring.delta_E[0, 0] / ring.delta_E[0, -1], 1,0.1), \
-        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[0, 0]} / {ring.delta_E[0, -1]}"
+        a = ring.delta_E[1,0]
+        b = ring.delta_E[1,1]
+        assert np.isclose(a/b,1,rtol=1e-1), \
+        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[1, 0]} == {ring.delta_E[1,1]}"
+
 
 
 
-- 
GitLab


From be1bf3f748f0f6c8aed50f4cb23508a90e2a5679 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Wed, 15 Jan 2025 13:31:34 +0100
Subject: [PATCH 25/48] re-added induced_voltage_mtw

---
 blond/impedances/impedance.py | 21 +++++++++++++++++++++
 1 file changed, 21 insertions(+)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 154aab7d..04772213 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -1044,6 +1044,9 @@ class InducedVoltageResonator(_InducedVoltage):
     def on_profile_change(self):
         self.process(self)
 
+
+
+
     def induced_voltage_1turn(self, beam_spectrum_dict={}):
         r"""
         Method to calculate the induced voltage through linearly
@@ -1067,6 +1070,24 @@ class InducedVoltageResonator(_InducedVoltage):
                                                 self.induced_voltage,
                                                 bm.precision.real_t))
 
+    def induced_voltage_mtw(self, beam_spectrum_dict={}):
+
+        """
+        Induced voltage method for InducedVoltageResonator.
+        mtw_memory is shifted by one turn, setting the final values in the array as 0.
+        @fbatsch
+        """
+
+        # shift the entries in array by 1 t_rev and set to 0
+        self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.array_length))
+        # remove one turn length of memory
+        self.mtw_memory = self.mtw_memory[self.array_length:]
+        # Induced voltage of the current turn calculation
+        self.induced_voltage_1turn(beam_spectrum_dict)
+        # Add induced voltage of the current turn to the memory from previous
+        self.mtw_memory[:int(self.n_time)] += self.induced_voltage
+        self.induced_voltage = self.mtw_memory[:self.n_time]
+
     def to_gpu(self, recursive=True):
         """
         Transfer all necessary arrays to the GPU
-- 
GitLab


From b4d30c212b9df3473f403e51cf69d814711eaadb Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Tue, 28 Jan 2025 11:32:36 +0100
Subject: [PATCH 26/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Leonard Sebastian Thiele <leonard.sebastian.thiele@cern.ch>
---
 blond/impedances/impedance.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 04772213..af7e77aa 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -944,7 +944,7 @@ class InducedVoltageResonator(_InducedVoltage):
         If left out, the induced voltage is calculated at the times of the
         line density
     atLineDensityTimes : boolean
-        flag indicating if the induced voltage has to be computed for timeArray
+        flag indicating if the induced voltage has to be computed for time_array
         or for the line density
     n_time : int
         length of tArray
-- 
GitLab


From 5f8bd5d64c6fac1d140c0b9d2517d31d830f5da9 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:29:04 +0000
Subject: [PATCH 27/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Leonard Sebastian Thiele <leonard.sebastian.thiele@cern.ch>
---
 blond/impedances/impedance.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index af7e77aa..41476a43 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -422,7 +422,6 @@ class _InducedVoltage:
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
 
-
         # Setting to zero to the last part to remove the contribution from the
         # front wake
         self.induced_voltage[self.n_induced_voltage -
-- 
GitLab


From a790f91d0ba84f6fd4987b6dfedbcb990aedd284 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:29:19 +0000
Subject: [PATCH 28/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Simon Albright <simon.albright@cern.ch>
---
 blond/impedances/impedance.py | 2 --
 1 file changed, 2 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 41476a43..38e6b840 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -460,8 +460,6 @@ class _InducedVoltage:
                                     self.time_mtw, self.mtw_memory,
                                     left=0, right=0)
 
-
-
     def _track(self):
         """
         Tracking method
-- 
GitLab


From c9de341465831253f8e253580e6861e06bbfb99d Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:29:58 +0000
Subject: [PATCH 29/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Simon Albright <simon.albright@cern.ch>
---
 blond/impedances/impedance.py | 3 ---
 1 file changed, 3 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 38e6b840..b6fe9e16 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -1041,9 +1041,6 @@ class InducedVoltageResonator(_InducedVoltage):
     def on_profile_change(self):
         self.process(self)
 
-
-
-
     def induced_voltage_1turn(self, beam_spectrum_dict={}):
         r"""
         Method to calculate the induced voltage through linearly
-- 
GitLab


From 9425dd9e6379336358b04d4693f176a68e41855e Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:32:00 +0000
Subject: [PATCH 30/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Simon Albright <simon.albright@cern.ch>
---
 unittests/input_parameters/test_ring.py | 3 ---
 1 file changed, 3 deletions(-)

diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index 197d79ae..facebc16 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -191,9 +191,6 @@ class TestGeneralParameters(unittest.TestCase):
                     [np.float64(1.0), np.float64(2.0), np.float64(3.0)])  # might crash
         Ring(C, alpha, momentum, Proton(), n_turns)  # might crash
 
-
-
-
     def test_bug_multi_RF_dE(self):
         from blond.beam.beam import Proton
         C = 2 * np.pi * 1100.009  # Ring circumference [m]
-- 
GitLab


From 886d90effdf3cf84db97fcdcc7f8c1589336551e Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:32:06 +0000
Subject: [PATCH 31/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Simon Albright <simon.albright@cern.ch>
---
 unittests/impedances/test_impedance.py | 21 ++++++++++++++-------
 1 file changed, 14 insertions(+), 7 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index f081af80..2212be9c 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -277,15 +277,18 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
 
 
     def test_init(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator)
 
     def test_init_mtw(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator,
                                       multi_turn_wake=True
                                       )
         
     def test_mtw_false_induced_volt(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator,
                                       multi_turn_wake=False)
         ivr.process()
         ivr.induced_voltage_1turn()
@@ -293,16 +296,20 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         self.assertTrue(np.any(my_array != 0.0))
 
     def test_multi_rf_station(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
-                                      rf_station=self.rf_station, multi_turn_wake=False)
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator, 
+                                      rf_station=self.rf_station,
+                                      multi_turn_wake=False)
         ivr.process()
         ivr.induced_voltage_1turn()
         my_array = ivr.induced_voltage
         self.assertTrue(np.any(my_array != 0.0))
 
     def test_mtw_true_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station, multi_turn_wake=True)
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator,
+                                      rf_station=self.rf_station,
+                                      multi_turn_wake=True)
         ivr.process()
         ivr.induced_voltage_mtw()
         my_array = ivr.induced_voltage
-- 
GitLab


From 629e67146f4081de63be4e2a32de7d8c62fb3450 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:32:27 +0000
Subject: [PATCH 32/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Simon Albright <simon.albright@cern.ch>
---
 unittests/impedances/test_impedance.py | 41 ++++++++++++--------------
 1 file changed, 19 insertions(+), 22 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 2212be9c..b882126e 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -237,28 +237,25 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         alpha_c = np.full(n_rf_stations, mom_compaction)
         n_turns = math.ceil(ring_per_section.n_turns /n_rf_stations)
 
-        ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c, particle=Proton(), n_turns=n_turns, n_sections=n_rf_stations,
-                    synchronous_data=25.92e9)
-
-        self.rf_station = RFStation(
-            ring=ring,
-            harmonic=[4620],
-            voltage=[0.9e6],
-            phi_rf_d=[0.0],
-            n_rf=1,
-        )
-        
-        beam = Beam(
-            ring=ring,
-            n_macroparticles=1001,
-            intensity=int(1e10),
-        )
-        bigaussian(
-            ring=ring,
-            rf_station=self.rf_station,
-            beam=beam,
-            sigma_dt=2e-9 / 4, seed=1,
-        )
+        ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c,
+                    particle=Proton(), n_turns=n_turns,
+                    n_sections=n_rf_stations, synchronous_data=25.92e9)
+
+        self.rf_station = RFStation(ring=ring,
+                                    harmonic=[4620],
+                                    voltage=[0.9e6],
+                                    phi_rf_d=[0.0],
+                                    n_rf=1)
+
+        beam = Beam(ring=ring,
+                    n_macroparticles=1001,
+                    intensity=int(1e10))
+
+        bigaussian(ring=ring,
+                   rf_station=self.rf_station,
+                   beam=beam,
+                   sigma_dt=2e-9 / 4,
+                   seed=1)
         self.beam = beam
         cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
                                  rf_station=self.rf_station, cuts_unit='rad')
-- 
GitLab


From 2736bd3782779a4980ebe0b0a486570aead51637 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:36:30 +0000
Subject: [PATCH 33/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Leonard Sebastian Thiele <leonard.sebastian.thiele@cern.ch>
---
 blond/impedances/impedance.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index b6fe9e16..c930dffc 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -427,7 +427,6 @@ class _InducedVoltage:
         self.induced_voltage[self.n_induced_voltage -
                              self.front_wake_buffer:] = 0
 
-
         # Add the induced voltage of the current turn to the memory from
         # previous turns
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
-- 
GitLab


From 91ef994a80d5563ecff7d934e17b80f16c26969d Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elleanor.rose.lamb@cern.ch>
Date: Thu, 6 Feb 2025 10:40:18 +0000
Subject: [PATCH 34/48] Apply 1 suggestion(s) to 1 file(s)

Co-authored-by: Leonard Sebastian Thiele <leonard.sebastian.thiele@cern.ch>
---
 unittests/impedances/test_impedance.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index b882126e..2e609315 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -280,8 +280,7 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
     def test_init_mtw(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
                                       resonators=self.resonator,
-                                      multi_turn_wake=True
-                                      )
+                                      multi_turn_wake=True)
         
     def test_mtw_false_induced_volt(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-- 
GitLab


From ae5c22599942b86ac111a4f5ae9a021a84f95776 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Tue, 4 Mar 2025 09:49:37 +0100
Subject: [PATCH 35/48] resolve merge request comments

---
 blond/impedances/impedance.py           |  24 ++---
 blond/input_parameters/ring.py          |  18 ++--
 unittests/impedances/test_impedance.py  | 120 +++++++++---------------
 unittests/input_parameters/test_ring.py |  16 ++--
 4 files changed, 70 insertions(+), 108 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index c930dffc..aa3ce3a6 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -412,9 +412,7 @@ class _InducedVoltage:
         from previous passages (multi-turn wake)
         """
 
-
         if beam_spectrum_dict is None:
-            print('beam_spectrum_dict is None')
             beam_spectrum_dict = dict()
         # Shift of the memory wake field by the current revolution period
         self.shift_trev()
@@ -431,6 +429,7 @@ class _InducedVoltage:
         # previous turns
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
+
     def shift_trev_freq(self):
         """
         Method to shift the induced voltage by a revolution period in the
@@ -935,7 +934,7 @@ class InducedVoltageResonator(_InducedVoltage):
         Copy of the Beam object in order to access the beam info.
     profile : Profile
         Copy of the Profile object in order to access the line density.
-    tArray : float array
+    timr_array : float array
         Array of time values where the induced voltage is calculated.
         If left out, the induced voltage is calculated at the times of the
         line density
@@ -943,7 +942,7 @@ class InducedVoltageResonator(_InducedVoltage):
         flag indicating if the induced voltage has to be computed for time_array
         or for the line density
     n_time : int
-        length of tArray
+        length of time_array
     R, omega_r, Q : lists of float
         Resonators parameters
     n_resonators : int
@@ -962,6 +961,7 @@ class InducedVoltageResonator(_InducedVoltage):
                  time_array: Optional[NDArray] = None,
                  rf_station: Optional[RFStation] = None,
                  array_length: Optional[int] = None,
+                 mtw_mode: Optional[MtwModeTypes] = None,
                  use_regular_fft: bool = True) -> None:
 
         # Test if one or more quality factors is smaller than 0.5.
@@ -978,15 +978,15 @@ class InducedVoltageResonator(_InducedVoltage):
         # If left out, the induced voltage is calculated at the times of the
         # line density.
         if time_array is None:
-            self.tArray = self.profile.bin_centers
+            self.time_array = self.profile.bin_centers
             self.atLineDensityTimes = True
         else:
-            self.tArray = time_array
+            self.time_array = time_array
             self.atLineDensityTimes = False
 
         self.array_length = array_length
 
-        self.n_time = len(self.tArray)
+        self.n_time = len(self.time_array)
 
         # Copy of the shunt impedances of the Resonators in* :math:`\Omega`
         self.R = resonators.R_S
@@ -1011,7 +1011,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self._kappa1 = np.zeros(
             int(self.profile.n_slices - 1), dtype=bm.precision.real_t, order='C')
 
-        # Matrix to hold n_times many tArray[t]-bin_centers arrays.
+        # Matrix to hold n_times many time_array[t]-bin_centers arrays.
         self._deltaT = np.zeros(
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
@@ -1020,7 +1020,7 @@ class InducedVoltageResonator(_InducedVoltage):
                          frequency_resolution=frequency_resolution,
                          wake_length=wake_length,
                          multi_turn_wake=multi_turn_wake,
-                         rf_station=rf_station, mtw_mode=None,
+                         rf_station=rf_station, mtw_mode=mtw_mode,
                          use_regular_fft=use_regular_fft)
 
     def process(self):
@@ -1052,7 +1052,7 @@ class InducedVoltageResonator(_InducedVoltage):
                                                 self.profile.bin_centers,
                                                 self.profile.bin_size,
                                                 self.n_time, self._deltaT,
-                                                self.tArray, self._reOmegaP,
+                                                self.time_array, self._reOmegaP,
                                                 self._imOmegaP, self._Qtilde,
                                                 self.n_resonators,
                                                 self.omega_r,
@@ -1093,7 +1093,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.array(self.induced_voltage)
         self._kappa1 = cp.array(self._kappa1)
         self._deltaT = cp.array(self._deltaT)
-        self.tArray = cp.array(self.tArray)
+        self.time_array = cp.array(self.time_array)
         self._tmp_matrix = cp.array(self._tmp_matrix)
         # to make sure it will not be called again
         self._device: DeviceType = 'GPU'
@@ -1110,7 +1110,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.asnumpy(self.induced_voltage)
         self._kappa1 = cp.asnumpy(self._kappa1)
         self._deltaT = cp.asnumpy(self._deltaT)
-        self.tArray = cp.asnumpy(self.tArray)
+        self.time_array = cp.asnumpy(self.time_array)
         self._tmp_matrix = cp.asnumpy(self._tmp_matrix)
 
         # to make sure it will not be called again
diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 7aae56c9..007f7998 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -276,27 +276,23 @@ class Ring:
         self.gamma: NDArray = np.sqrt(1 + (self.momentum / self.particle.mass) ** 2)
         self.energy: NDArray = np.sqrt(self.momentum ** 2 + self.particle.mass ** 2)
         self.kin_energy: NDArray = (np.sqrt(self.momentum ** 2 + self.particle.mass ** 2) -
-                                       self.particle.mass)
+                                    self.particle.mass)
         self.t_rev: NDArray = np.dot(self.ring_length, 1 / (self.beta * c))
         self.cycle_time: NDArray = np.cumsum(self.t_rev)  # Always starts with zero
         self.f_rev: NDArray = 1 / self.t_rev
         self.omega_rev: NDArray = 2 * np.pi * self.f_rev
 
-        
-
         if self.n_sections == 1:
             self.delta_E = np.diff(self.energy, axis=1)
-
         else:
-            self.delta_E = np.reshape(np.zeros(n_turns*n_sections),(n_sections,n_turns))
-            for i in range(n_turns):  
-                for j in range(n_sections):
+            # when there is more than 1 RF station, self.energy has shape (n_sections, n_turns+1)
+            self.delta_E = np.zeros((n_sections, n_turns))
+            for j in range(n_sections):
+                for i in range(n_turns):
                     if j == 0:
-                        self.delta_E[j,i] = self.energy[j,i+1]-self.energy[-1,i]
+                        self.delta_E[j, i] = self.energy[j, i + 1] - self.energy[-1, i]
                     else:
-                        self.delta_E[j,i] = self.energy[j,i+1]-self.energy[j-1,i+1]
-
-
+                        self.delta_E[j, i] = self.energy[j, i + 1] - self.energy[j - 1, i + 1]
 
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 2e609315..e8926293 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -129,61 +129,54 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
         from blond.impedances.impedance_sources import Resonators
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
-        import math
 
-        gamma_transition = 1 / np.sqrt(0.00192)  
+        gamma_transition = 1 / np.sqrt(0.00192)
         mom_compaction = 1 / gamma_transition ** 2
         circ = 6911.56
 
-        ring = Ring(
-            ring_length=circ,
-            alpha_0=mom_compaction,
-            synchronous_data=25.92e9,
-            particle=Proton(),
-            n_turns=2)
-        
-        self.rf_station = RFStation(
-            ring=ring,
-            harmonic=[4620],
-            voltage=[0.9e6],
-            phi_rf_d=[0.0],
-            n_rf=1,
-        )
-
-        beam = Beam(
-            ring=ring,
-            n_macroparticles=1001,
-            intensity=int(1e10),
-        )
-        bigaussian(
-            ring=ring,
-            rf_station=self.rf_station,
-            beam=beam,
-            sigma_dt=2e-9 / 4, seed=1,
-        )
+        ring = Ring(ring_length=circ,
+                    alpha_0=mom_compaction,
+                    synchronous_data=25.92e9,
+                    particle=Proton(),
+                    n_turns=2)
+
+        self.rf_station = RFStation(ring=ring,
+                                    harmonic=[4620],
+                                    voltage=[0.9e6],
+                                    phi_rf_d=[0.0],
+                                    n_rf=1,
+                                    )
+
+        beam = Beam(ring=ring,
+                    n_macroparticles=1001,
+                    intensity=int(1e10),
+                    )
+        bigaussian(ring=ring,
+                   rf_station=self.rf_station,
+                   beam=beam,
+                   sigma_dt=2e-9 / 4, seed=1,
+                   )
         self.beam = beam
 
         cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
                                  rf_station=self.rf_station, cuts_unit='rad')
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        
+
         this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10**6
-        self.f_res = np.load(this_directory + './data/f_res.npy') * 10**9
+        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10 ** 6
+        self.f_res = np.load(this_directory + './data/f_res.npy') * 10 ** 9
         self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
 
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
-
     def test_init(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
 
     def test_init_mtw(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=True
-                                      )
-        
+                                      multi_turn_wake=True)
+
     def test_mtw_false_induced_volt(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       multi_turn_wake=False)
@@ -193,7 +186,7 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
         self.assertTrue(np.any(my_array != 0.0))
 
     def test_multi_rf_station(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator, 
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       rf_station=self.rf_station, multi_turn_wake=False)
         ivr.process()
         ivr.induced_voltage_1turn()
@@ -209,70 +202,52 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
         self.assertTrue(np.any(my_array != 0.0))
 
 
-
-
 class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
 
     def setUp(self):
         import os
         import numpy as np
         from blond.beam.beam import Beam, Proton
-        from blond.beam.distributions import bigaussian
         from blond.beam.profile import CutOptions, FitOptions, Profile
         from blond.impedances.impedance_sources import Resonators
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
-        import math 
-        
-        gamma_transition = 1 / np.sqrt(0.00192)  
+
+        gamma_transition = 1 / np.sqrt(0.00192)
         mom_compaction = 1 / gamma_transition ** 2
         circ = 6911.56
         n_rf_stations = 2
-        l_per_section = circ/ n_rf_stations
-        ring_per_section = Ring(ring_length=l_per_section, 
-                                alpha_0=mom_compaction,
-                        particle=Proton(), synchronous_data=25.92e9)
-
+        l_per_section = circ / n_rf_stations
+        n_turns = 2
         sectionlengths = np.full(n_rf_stations, l_per_section)
         alpha_c = np.full(n_rf_stations, mom_compaction)
-        n_turns = math.ceil(ring_per_section.n_turns /n_rf_stations)
+        energy_program = np.full((n_rf_stations, n_turns + 1), 25.92e9)
 
         ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c,
                     particle=Proton(), n_turns=n_turns,
-                    n_sections=n_rf_stations, synchronous_data=25.92e9)
-
+                    n_sections=n_rf_stations, synchronous_data=energy_program)
         self.rf_station = RFStation(ring=ring,
-                                    harmonic=[4620],
-                                    voltage=[0.9e6],
-                                    phi_rf_d=[0.0],
-                                    n_rf=1)
+                                     harmonic=[4620],
+                                     voltage=[0.9e6],
+                                     phi_rf_d=[0.0],
+                                     n_rf=1)
 
         beam = Beam(ring=ring,
                     n_macroparticles=1001,
                     intensity=int(1e10))
 
-        bigaussian(ring=ring,
-                   rf_station=self.rf_station,
-                   beam=beam,
-                   sigma_dt=2e-9 / 4,
-                   seed=1)
         self.beam = beam
         cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
-                                 rf_station=self.rf_station, cuts_unit='rad')
+                                 cuts_unit='rad', rf_station=self.rf_station)
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
-        
+
         this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10**6
-        self.f_res = np.load(this_directory + './data/f_res.npy') * 10**9
+        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10 ** 6
+        self.f_res = np.load(this_directory + './data/f_res.npy') * 10 ** 9
         self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
-
-
-
-
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
-
     def test_init(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
                                       resonators=self.resonator)
@@ -281,7 +256,7 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
                                       resonators=self.resonator,
                                       multi_turn_wake=True)
-        
+
     def test_mtw_false_induced_volt(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
                                       resonators=self.resonator,
@@ -293,7 +268,7 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
 
     def test_multi_rf_station(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator, 
+                                      resonators=self.resonator,
                                       rf_station=self.rf_station,
                                       multi_turn_wake=False)
         ivr.process()
@@ -312,10 +287,5 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         self.assertTrue(np.any(my_array != 0.0))
 
 
-
-
-
-
-
 if __name__ == '__main__':
     unittest.main()
diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index facebc16..11834571 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -198,20 +198,16 @@ class TestGeneralParameters(unittest.TestCase):
         alpha = 1 / gamma_t ** 2  # Momentum compaction factor
         n_turns = self.n_turns
         n_sections = 4
-        l_per_section = C/ n_sections
+        l_per_section = C / n_sections
         section_lengths = np.full(n_sections, l_per_section)
         # should be a linear energy ramp, taking the value from the section in the previous turn
-        momentum = np.arange(n_sections * (n_turns + 1),dtype=float).reshape(n_turns + 1, n_sections).T*1e10
-        momentum[:, 0] = momentum[0, 1]
+        momentum = np.arange(0, (n_turns+1)*n_sections).reshape((n_turns+1, n_sections)).transpose()+1e10
         ring = Ring(ring_length=section_lengths, alpha_0=alpha,
-            synchronous_data=momentum, particle=Proton(),
-            n_turns=n_turns, n_sections=n_sections,
+                    synchronous_data=momentum, particle=Proton(),
+                    n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
-        a = ring.delta_E[1,0]
-        b = ring.delta_E[1,1]
-        assert np.isclose(a/b,1,rtol=1e-1), \
-        f"Assertion failed for linear momentum with multi RF: {ring.delta_E[1, 0]} == {ring.delta_E[1,1]}"
-
+        delta_E_values = ring.delta_E.flatten()
+        np.testing.assert_allclose(delta_E_values, delta_E_values[0], rtol=1e-2,)
 
 
 
-- 
GitLab


From 66917d4242223ae7869a9b21f342215c030d66c6 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Tue, 4 Mar 2025 13:10:35 +0100
Subject: [PATCH 36/48] resolve merge request comments

---
 unittests/input_parameters/test_ring.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index 11834571..fd312c8c 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -201,7 +201,7 @@ class TestGeneralParameters(unittest.TestCase):
         l_per_section = C / n_sections
         section_lengths = np.full(n_sections, l_per_section)
         # should be a linear energy ramp, taking the value from the section in the previous turn
-        momentum = np.arange(0, (n_turns+1)*n_sections).reshape((n_turns+1, n_sections)).transpose()+1e10
+        momentum = np.arange(0, (n_turns+1)*n_sections).reshape((n_turns+1, n_sections)).transpose()*1e10
         ring = Ring(ring_length=section_lengths, alpha_0=alpha,
                     synchronous_data=momentum, particle=Proton(),
                     n_turns=n_turns, n_sections=n_sections,
-- 
GitLab


From ca73992674c1ddaf1097803550c6ce7eaee981a7 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Tue, 4 Mar 2025 13:13:46 +0100
Subject: [PATCH 37/48] resolve merge request comments

---
 unittests/input_parameters/test_ring.py | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index fd312c8c..47f26b40 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -22,6 +22,7 @@ from blond.input_parameters.ring import Ring
 from blond.input_parameters.ring_options import convert_data
 
 
+
 class TestGeneralParameters(unittest.TestCase):
 
     # Initialization ----------------------------------------------------------
@@ -193,6 +194,7 @@ class TestGeneralParameters(unittest.TestCase):
 
     def test_bug_multi_RF_dE(self):
         from blond.beam.beam import Proton
+        # test if delta_E is calculated correctly for multi-RF stations
         C = 2 * np.pi * 1100.009  # Ring circumference [m]
         gamma_t = 18.0  # Gamma at transition
         alpha = 1 / gamma_t ** 2  # Momentum compaction factor
@@ -200,16 +202,14 @@ class TestGeneralParameters(unittest.TestCase):
         n_sections = 4
         l_per_section = C / n_sections
         section_lengths = np.full(n_sections, l_per_section)
-        # should be a linear energy ramp, taking the value from the section in the previous turn
-        momentum = np.arange(0, (n_turns+1)*n_sections).reshape((n_turns+1, n_sections)).transpose()*1e10
+        # a linear energy ramp, taking the value from the section in the previous turn
+        momentum = np.arange(0, (n_turns + 1) * n_sections).reshape((n_turns + 1, n_sections)).transpose() * 1e10
         ring = Ring(ring_length=section_lengths, alpha_0=alpha,
                     synchronous_data=momentum, particle=Proton(),
                     n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
         delta_E_values = ring.delta_E.flatten()
-        np.testing.assert_allclose(delta_E_values, delta_E_values[0], rtol=1e-2,)
-
-
+        np.testing.assert_allclose(delta_E_values, delta_E_values[0], rtol=1e-2)
 
 
 if __name__ == '__main__':
-- 
GitLab


From 80106d4371bafde960bbe2522dad7815b5a8b148 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Tue, 4 Mar 2025 13:18:58 +0100
Subject: [PATCH 38/48] resolve merge request comments

---
 unittests/input_parameters/test_ring.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index 47f26b40..2665906a 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -209,7 +209,7 @@ class TestGeneralParameters(unittest.TestCase):
                     n_turns=n_turns, n_sections=n_sections,
                     synchronous_data_type='momentum')
         delta_E_values = ring.delta_E.flatten()
-        np.testing.assert_allclose(delta_E_values, delta_E_values[0], rtol=1e-2)
+        np.testing.assert_allclose(delta_E_values, delta_E_values[0], rtol=1e-3)
 
 
 if __name__ == '__main__':
-- 
GitLab


From ec94be81d78e2d8c6f637ca4b2a45bb1df66fab6 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 7 Mar 2025 16:10:17 +0100
Subject: [PATCH 39/48] changed formatting

---
 blond/impedances/impedance.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index aa3ce3a6..389abc59 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -954,7 +954,7 @@ class InducedVoltageResonator(_InducedVoltage):
     @handle_legacy_kwargs
     def __init__(self, beam: Beam,
                  profile: Profile,
-                 resonators: Resonators,
+                 resonators: Optional[Resonators] = None,
                  frequency_resolution: Optional[float] = None,
                  wake_length: Optional[float] = None,
                  multi_turn_wake: bool = False,
-- 
GitLab


From e141141e1de83108106fc9cc632b58256178b7a6 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 7 Mar 2025 16:26:45 +0100
Subject: [PATCH 40/48] addded back time array as necessary to maintain length
 of mtw memory

---
 blond/impedances/impedance.py | 15 ++++++---------
 1 file changed, 6 insertions(+), 9 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 389abc59..c98785af 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -978,15 +978,15 @@ class InducedVoltageResonator(_InducedVoltage):
         # If left out, the induced voltage is calculated at the times of the
         # line density.
         if time_array is None:
-            self.time_array = self.profile.bin_centers
+            self.tArray = self.profile.bin_centers
             self.atLineDensityTimes = True
         else:
-            self.time_array = time_array
+            self.tArray = time_array
             self.atLineDensityTimes = False
 
         self.array_length = array_length
 
-        self.n_time = len(self.time_array)
+        self.n_time = len(self.tArray)
 
         # Copy of the shunt impedances of the Resonators in* :math:`\Omega`
         self.R = resonators.R_S
@@ -1037,9 +1037,6 @@ class InducedVoltageResonator(_InducedVoltage):
         self._deltaT = np.zeros(
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
-    def on_profile_change(self):
-        self.process(self)
-
     def induced_voltage_1turn(self, beam_spectrum_dict={}):
         r"""
         Method to calculate the induced voltage through linearly
@@ -1052,7 +1049,7 @@ class InducedVoltageResonator(_InducedVoltage):
                                                 self.profile.bin_centers,
                                                 self.profile.bin_size,
                                                 self.n_time, self._deltaT,
-                                                self.time_array, self._reOmegaP,
+                                                self.tArray, self._reOmegaP,
                                                 self._imOmegaP, self._Qtilde,
                                                 self.n_resonators,
                                                 self.omega_r,
@@ -1093,7 +1090,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.array(self.induced_voltage)
         self._kappa1 = cp.array(self._kappa1)
         self._deltaT = cp.array(self._deltaT)
-        self.time_array = cp.array(self.time_array)
+        self.tArray = cp.array(self.tArray)
         self._tmp_matrix = cp.array(self._tmp_matrix)
         # to make sure it will not be called again
         self._device: DeviceType = 'GPU'
@@ -1110,7 +1107,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.asnumpy(self.induced_voltage)
         self._kappa1 = cp.asnumpy(self._kappa1)
         self._deltaT = cp.asnumpy(self._deltaT)
-        self.time_array = cp.asnumpy(self.time_array)
+        self.tArray = cp.asnumpy(self.tArray)
         self._tmp_matrix = cp.asnumpy(self._tmp_matrix)
 
         # to make sure it will not be called again
-- 
GitLab


From 48c9c82fcbb3a7fff1988d08a77f88271b9678ac Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 7 Mar 2025 16:50:52 +0100
Subject: [PATCH 41/48] addded back time array as necessary to maintain length
 of mtw memory

---
 blond/impedances/impedance.py           | 10 ++++++++
 unittests/impedances/test_impedance.py  | 34 +++++++++++++++----------
 unittests/input_parameters/test_ring.py |  1 +
 3 files changed, 31 insertions(+), 14 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index c98785af..b7b4387d 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -424,11 +424,15 @@ class _InducedVoltage:
         # front wake
         self.induced_voltage[self.n_induced_voltage -
                              self.front_wake_buffer:] = 0
+        print(self.induced_voltage)
 
         # Add the induced voltage of the current turn to the memory from
         # previous turns
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
+        print(self.mtw_memory)
+
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
+        print(self.induced_voltage)
 
     def shift_trev_freq(self):
         """
@@ -1070,13 +1074,19 @@ class InducedVoltageResonator(_InducedVoltage):
 
         # shift the entries in array by 1 t_rev and set to 0
         self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.array_length))
+        print(self.mtw_memory)
         # remove one turn length of memory
         self.mtw_memory = self.mtw_memory[self.array_length:]
+        print(self.mtw_memory)
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
+        print(self.induced_voltage)
         # Add induced voltage of the current turn to the memory from previous
         self.mtw_memory[:int(self.n_time)] += self.induced_voltage
+        print(self.mtw_memory)
+
         self.induced_voltage = self.mtw_memory[:self.n_time]
+        print(self.induced_voltage)
 
     def to_gpu(self, recursive=True):
         """
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index e8926293..2fb4aaae 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -19,6 +19,7 @@ import numpy as np
 from blond.beam.profile import CutOptions, Profile
 from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime, InducedVoltageResonator
 from blond.impedances.impedance_sources import Resonators
+from impedances.impedance import TotalInducedVoltage
 
 
 class TestInducedVoltageFreq(unittest.TestCase):
@@ -190,16 +191,14 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
                                       rf_station=self.rf_station, multi_turn_wake=False)
         ivr.process()
         ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
+        self.assertTrue(np.any(ivr.induced_voltage != 0.0))
 
     def test_mtw_true_induced_voltage(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       rf_station=self.rf_station, multi_turn_wake=True)
         ivr.process()
         ivr.induced_voltage_mtw()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
+        self.assertTrue(np.any(ivr.induced_voltage != 0.0))
 
 
 class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
@@ -237,7 +236,7 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
                     intensity=int(1e10))
 
         self.beam = beam
-        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
+        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 9,
                                  cuts_unit='rad', rf_station=self.rf_station)
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
@@ -248,6 +247,21 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
 
+    def test_mtw_true_induced_voltage(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
+                                      resonators=self.resonator,
+                                      rf_station=self.rf_station,
+                                      multi_turn_wake=True)
+
+        #total_induced_voltage = TotalInducedVoltage(beam=beam, profile=profile, induced_voltage_list=[ivr])
+        #total_induced_voltage.sum()
+
+        #ivr.process()
+        #ivr.induced_voltage_mtw()
+        #my_array = ivr.induced_voltage
+        #self.assertTrue(np.any(my_array != 0.0))
+
+
     def test_init(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
                                       resonators=self.resonator)
@@ -276,15 +290,7 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         my_array = ivr.induced_voltage
         self.assertTrue(np.any(my_array != 0.0))
 
-    def test_mtw_true_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator,
-                                      rf_station=self.rf_station,
-                                      multi_turn_wake=True)
-        ivr.process()
-        ivr.induced_voltage_mtw()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
+
 
 
 if __name__ == '__main__':
diff --git a/unittests/input_parameters/test_ring.py b/unittests/input_parameters/test_ring.py
index 11834571..46918200 100644
--- a/unittests/input_parameters/test_ring.py
+++ b/unittests/input_parameters/test_ring.py
@@ -192,6 +192,7 @@ class TestGeneralParameters(unittest.TestCase):
         Ring(C, alpha, momentum, Proton(), n_turns)  # might crash
 
     def test_bug_multi_RF_dE(self):
+        # test if the delta_E is correctly calculated with multi-RF stations
         from blond.beam.beam import Proton
         C = 2 * np.pi * 1100.009  # Ring circumference [m]
         gamma_t = 18.0  # Gamma at transition
-- 
GitLab


From b23ceeee2d2ea949abff96d285d33529063237ee Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 7 Mar 2025 17:21:33 +0100
Subject: [PATCH 42/48] addded back time array as necessary to maintain length
 of mtw memory

---
 blond/impedances/impedance.py          |  3 +-
 unittests/impedances/test_impedance.py | 64 +++++++++++++-------------
 2 files changed, 34 insertions(+), 33 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index b7b4387d..7ab0e035 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -416,6 +416,7 @@ class _InducedVoltage:
             beam_spectrum_dict = dict()
         # Shift of the memory wake field by the current revolution period
         self.shift_trev()
+        print(self.induced_voltage)
 
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
@@ -460,7 +461,7 @@ class _InducedVoltage:
         # self.mtw_memory = bm.interp_const_space(self.time_mtw + t_rev,
         self.mtw_memory = bm.interp(self.time_mtw + t_rev,
                                     self.time_mtw, self.mtw_memory,
-                                    left=0, right=0)
+                                    left=0, right=0) # todo
 
     def _track(self):
         """
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 2fb4aaae..9acb2dae 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -19,7 +19,7 @@ import numpy as np
 from blond.beam.profile import CutOptions, Profile
 from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime, InducedVoltageResonator
 from blond.impedances.impedance_sources import Resonators
-from impedances.impedance import TotalInducedVoltage
+from blond.impedances.impedance import TotalInducedVoltage
 
 
 class TestInducedVoltageFreq(unittest.TestCase):
@@ -143,7 +143,7 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
 
         self.rf_station = RFStation(ring=ring,
                                     harmonic=[4620],
-                                    voltage=[0.9e6],
+                                    voltage=[3000e6],
                                     phi_rf_d=[0.0],
                                     n_rf=1,
                                     )
@@ -186,9 +186,18 @@ class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
         my_array = ivr.induced_voltage
         self.assertTrue(np.any(my_array != 0.0))
 
-    def test_multi_rf_station(self):
+    def test_total_induced_voltage(self):
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
                                       rf_station=self.rf_station, multi_turn_wake=False)
+        total_induced_voltage = TotalInducedVoltage(beam=self.beam, profile=self.profile, induced_voltage_list=[ivr])
+        total_induced_voltage.induced_voltage_sum()
+
+        #ivr.process()
+        #ivr.induced_voltage_mtw()
+        #my_array = ivr.induced_voltage
+        #self.assertTrue(np.any(my_array != 0.0))
+
+
         ivr.process()
         ivr.induced_voltage_1turn()
         self.assertTrue(np.any(ivr.induced_voltage != 0.0))
@@ -220,14 +229,15 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         n_turns = 2
         sectionlengths = np.full(n_rf_stations, l_per_section)
         alpha_c = np.full(n_rf_stations, mom_compaction)
-        energy_program = np.full((n_rf_stations, n_turns + 1), 25.92e9)
+        energy_program = np.array([[1*25e9,1*25e9,1.2*25e9],
+                                   [1.0*25e9,1.1*25e9,1.3*25e9]]) # todo
 
         ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c,
                     particle=Proton(), n_turns=n_turns,
                     n_sections=n_rf_stations, synchronous_data=energy_program)
         self.rf_station = RFStation(ring=ring,
                                      harmonic=[4620],
-                                     voltage=[0.9e6],
+                                     voltage=[300e6],
                                      phi_rf_d=[0.0],
                                      n_rf=1)
 
@@ -241,25 +251,12 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         self.profile = Profile(beam, cut_options,
                                FitOptions(fit_option='gaussian'))
 
-        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10 ** 6
-        self.f_res = np.load(this_directory + './data/f_res.npy') * 10 ** 9
-        self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
+        self.R_shunt = 11897424000
+        self.f_res = 11897424000.
+        self.Q_factor = 696000.0
         self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
-
-    def test_mtw_true_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator,
-                                      rf_station=self.rf_station,
-                                      multi_turn_wake=True)
-
-        #total_induced_voltage = TotalInducedVoltage(beam=beam, profile=profile, induced_voltage_list=[ivr])
-        #total_induced_voltage.sum()
-
-        #ivr.process()
-        #ivr.induced_voltage_mtw()
-        #my_array = ivr.induced_voltage
-        #self.assertTrue(np.any(my_array != 0.0))
+        print(self.Q_factor)
+        print(self.f_res)
 
 
     def test_init(self):
@@ -280,15 +277,18 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         my_array = ivr.induced_voltage
         self.assertTrue(np.any(my_array != 0.0))
 
-    def test_multi_rf_station(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator,
-                                      rf_station=self.rf_station,
-                                      multi_turn_wake=False)
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
+
+    def test_total_induced_voltage(self):
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
+                                      rf_station=self.rf_station, multi_turn_wake=False,
+                                      wake_length=self.profile.bin_size * len(self.profile.bin_centers),
+                                      array_length=2**9)
+
+        total_induced_voltage = TotalInducedVoltage(beam=self.beam, profile=self.profile,
+                                                    induced_voltage_list=[ivr], )
+        total_induced_voltage.induced_voltage_sum() # todo
+        total_induced_voltage.induced_voltage_sum()
+
 
 
 
-- 
GitLab


From 055b8b4bd6e4bb52f8589fef92f42cbbaf2b8ff1 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Thu, 13 Mar 2025 09:58:49 +0100
Subject: [PATCH 43/48] merge request comments

---
 unittests/impedances/data/Q_factor.npy | Bin 544 -> 0 bytes
 unittests/impedances/data/R_shunt.npy  | Bin 544 -> 0 bytes
 unittests/impedances/data/f_res.npy    | Bin 544 -> 0 bytes
 3 files changed, 0 insertions(+), 0 deletions(-)
 delete mode 100644 unittests/impedances/data/Q_factor.npy
 delete mode 100644 unittests/impedances/data/R_shunt.npy
 delete mode 100644 unittests/impedances/data/f_res.npy

diff --git a/unittests/impedances/data/Q_factor.npy b/unittests/impedances/data/Q_factor.npy
deleted file mode 100644
index 1132114819c01ac1a81b9235d42d363a97760e18..0000000000000000000000000000000000000000
GIT binary patch
literal 0
HcmV?d00001

literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<&-aHw|x(+nyQ8pemw3@)hR40=%Y6;Sa~DBlD^E7U;w
z4G{W36NJ|2gwPHXp!_}v&EN^41ExS|n0kgx2>(GVDh)FqrVd6Iv_sT6^h0TwzwA)O
zVgB4u3=xOfb07%9e^3CS(fKfS1`{Fj4c-uX0n8tM5PksEpAS+Y{27xV^ou=Ux*`7n
Ygs-z7LO1P$&@8(k^t){kx@VsQ02D)1r2qf`

diff --git a/unittests/impedances/data/R_shunt.npy b/unittests/impedances/data/R_shunt.npy
deleted file mode 100644
index f45cbe077aab4a0896cd1cd41bde6b2ca7ccbd04..0000000000000000000000000000000000000000
GIT binary patch
literal 0
HcmV?d00001

literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<(2Xjmui2%;GdG&z8%2Ck1_n!)QUgnn}wOgGH=1>qk%
z2%#0`Kxm)65ISH9gpRlkp@Z%~<^O`|2G7qBe!+97{22({u^B=)9f8nUC!peoA@u(z
z5c=m+2wipzLQmNUp}*{c&?4s`bj?u+9kU5S&pH65FF<IheF1^WV1GX-kpt5W=3-#l
zVKyt6X88RUqTcu}n0}xL38w=qc_H%uAm$zTs0iVogs3lQR|oSCTxJ8)3D-5i^a3?O
cFx?Q%0Hzt5p#GZm0b)<eO9&0~C+iaj03BFuSpWb4

diff --git a/unittests/impedances/data/f_res.npy b/unittests/impedances/data/f_res.npy
deleted file mode 100644
index 72fff53c63ff2c43e265be3a912da52b2d4479b2..0000000000000000000000000000000000000000
GIT binary patch
literal 0
HcmV?d00001

literal 544
zcmbR27wQ`j$;eQ~P_3SlTAW;@Zl$1ZlV+i=qoAIaUsO_*m=~X4l#&V(cT3DEP6dh=
zXCxM+0{I%IMmm~03bhL411<&zh6RuJ9dZQH3>!I4f@p{LN6v%kxf~ZjbOXn?>md5T
zW1%}B`a*r%0}!1NaQ_jAj$ru(rVH$LJpu6<{4$?`=!zuf7a)4U?s72wL-p)S5PwF6
z%4-mPqLcLvh<>nw?=6UCNIeUt8-Az01Mvku#J&g74p4U-&`JTjXTxsAyC89c`q%eB
zbV7UoeGpyYlmvEHLx|yH5PyekJ=i@L4oidWWocXe6(kO$7o1?f1mZJX5V{PaKb#i5
zijUsVaPkUBULfTGm{t%@1&7Omm4Css!+GT!Ao&YwyKjSNh2xsw@Y->88aR9y`1c(F
mi67W6cN#=5xO5vF&Ixt`5P$tNy$BL#*cuIXSAjC~6-NMe6rXkg

-- 
GitLab


From 6628534efec406d854bd2ca1d6f41e14a12da2d5 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Thu, 13 Mar 2025 09:58:59 +0100
Subject: [PATCH 44/48] merge request comments

---
 blond/impedances/impedance.py          |  23 +--
 unittests/impedances/test_impedance.py | 185 ++++++-------------------
 2 files changed, 45 insertions(+), 163 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 7ab0e035..6c280980 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -411,12 +411,10 @@ class _InducedVoltage:
         Method to calculate the induced voltage taking into account the effect
         from previous passages (multi-turn wake)
         """
-
         if beam_spectrum_dict is None:
             beam_spectrum_dict = dict()
         # Shift of the memory wake field by the current revolution period
         self.shift_trev()
-        print(self.induced_voltage)
 
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
@@ -425,15 +423,12 @@ class _InducedVoltage:
         # front wake
         self.induced_voltage[self.n_induced_voltage -
                              self.front_wake_buffer:] = 0
-        print(self.induced_voltage)
 
         # Add the induced voltage of the current turn to the memory from
         # previous turns
         self.mtw_memory[:self.n_induced_voltage] += self.induced_voltage
-        print(self.mtw_memory)
 
         self.induced_voltage = self.mtw_memory[:self.n_induced_voltage]
-        print(self.induced_voltage)
 
     def shift_trev_freq(self):
         """
@@ -457,11 +452,13 @@ class _InducedVoltage:
         """
 
         t_rev = self.rf_params.t_rev[self.rf_params.counter[0]]
+        print(t_rev)
+
 
         # self.mtw_memory = bm.interp_const_space(self.time_mtw + t_rev,
         self.mtw_memory = bm.interp(self.time_mtw + t_rev,
                                     self.time_mtw, self.mtw_memory,
-                                    left=0, right=0) # todo
+                                    left=0, right=0)  # todo
 
     def _track(self):
         """
@@ -1001,7 +998,6 @@ class InducedVoltageResonator(_InducedVoltage):
         self.Q = resonators.Q
         # Number of resonators
         self.n_resonators = len(self.R)
-
         # For internal use
         self._Qtilde = self.Q * np.sqrt(1. - 1. / (4. * self.Q ** 2.))
         self._reOmegaP = self.omega_r * self._Qtilde / self.Q
@@ -1064,7 +1060,7 @@ class InducedVoltageResonator(_InducedVoltage):
                                                 self.beam.ratio, self.R,
                                                 self.induced_voltage,
                                                 bm.precision.real_t))
-
+    '''
     def induced_voltage_mtw(self, beam_spectrum_dict={}):
 
         """
@@ -1072,23 +1068,16 @@ class InducedVoltageResonator(_InducedVoltage):
         mtw_memory is shifted by one turn, setting the final values in the array as 0.
         @fbatsch
         """
-
         # shift the entries in array by 1 t_rev and set to 0
         self.mtw_memory = np.append(self.mtw_memory, np.zeros(self.array_length))
-        print(self.mtw_memory)
         # remove one turn length of memory
         self.mtw_memory = self.mtw_memory[self.array_length:]
-        print(self.mtw_memory)
         # Induced voltage of the current turn calculation
         self.induced_voltage_1turn(beam_spectrum_dict)
-        print(self.induced_voltage)
         # Add induced voltage of the current turn to the memory from previous
         self.mtw_memory[:int(self.n_time)] += self.induced_voltage
-        print(self.mtw_memory)
-
         self.induced_voltage = self.mtw_memory[:self.n_time]
-        print(self.induced_voltage)
-
+    '''
     def to_gpu(self, recursive=True):
         """
         Transfer all necessary arrays to the GPU
@@ -1122,4 +1111,4 @@ class InducedVoltageResonator(_InducedVoltage):
         self._tmp_matrix = cp.asnumpy(self._tmp_matrix)
 
         # to make sure it will not be called again
-        self._device: DeviceType = 'CPU'
+        self._device: DeviceType = 'CPU'
\ No newline at end of file
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 9acb2dae..e96ef512 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -119,98 +119,7 @@ class TestInducedVoltageTime(unittest.TestCase):
         np.testing.assert_allclose(test_object.wake_length_input, 11e-9)
 
 
-class TestInducedVoltageResonator_n_rf1(unittest.TestCase):
-
-    def setUp(self):
-        import os
-        import numpy as np
-        from blond.beam.beam import Beam, Proton
-        from blond.beam.distributions import bigaussian
-        from blond.beam.profile import CutOptions, FitOptions, Profile
-        from blond.impedances.impedance_sources import Resonators
-        from blond.input_parameters.rf_parameters import RFStation
-        from blond.input_parameters.ring import Ring
-
-        gamma_transition = 1 / np.sqrt(0.00192)
-        mom_compaction = 1 / gamma_transition ** 2
-        circ = 6911.56
-
-        ring = Ring(ring_length=circ,
-                    alpha_0=mom_compaction,
-                    synchronous_data=25.92e9,
-                    particle=Proton(),
-                    n_turns=2)
-
-        self.rf_station = RFStation(ring=ring,
-                                    harmonic=[4620],
-                                    voltage=[3000e6],
-                                    phi_rf_d=[0.0],
-                                    n_rf=1,
-                                    )
-
-        beam = Beam(ring=ring,
-                    n_macroparticles=1001,
-                    intensity=int(1e10),
-                    )
-        bigaussian(ring=ring,
-                   rf_station=self.rf_station,
-                   beam=beam,
-                   sigma_dt=2e-9 / 4, seed=1,
-                   )
-        self.beam = beam
-
-        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 8,
-                                 rf_station=self.rf_station, cuts_unit='rad')
-        self.profile = Profile(beam, cut_options,
-                               FitOptions(fit_option='gaussian'))
-
-        this_directory = os.path.dirname(os.path.realpath(__file__)) + '/'
-        self.R_shunt = np.load(this_directory + './data/R_shunt.npy') * 10 ** 6
-        self.f_res = np.load(this_directory + './data/f_res.npy') * 10 ** 9
-        self.Q_factor = np.load(this_directory + './data/Q_factor.npy')
-
-        self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
-
-    def test_init(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator)
-
-    def test_init_mtw(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=True)
-
-    def test_mtw_false_induced_volt(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      multi_turn_wake=False)
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
-
-    def test_total_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station, multi_turn_wake=False)
-        total_induced_voltage = TotalInducedVoltage(beam=self.beam, profile=self.profile, induced_voltage_list=[ivr])
-        total_induced_voltage.induced_voltage_sum()
-
-        #ivr.process()
-        #ivr.induced_voltage_mtw()
-        #my_array = ivr.induced_voltage
-        #self.assertTrue(np.any(my_array != 0.0))
-
-
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        self.assertTrue(np.any(ivr.induced_voltage != 0.0))
-
-    def test_mtw_true_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station, multi_turn_wake=True)
-        ivr.process()
-        ivr.induced_voltage_mtw()
-        self.assertTrue(np.any(ivr.induced_voltage != 0.0))
-
-
-class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
+class TestInducedVoltageResonator_multi_turn_wake(unittest.TestCase):
 
     def setUp(self):
         import os
@@ -220,78 +129,62 @@ class TestInducedVoltageResonator_nrf_2(unittest.TestCase):
         from blond.impedances.impedance_sources import Resonators
         from blond.input_parameters.rf_parameters import RFStation
         from blond.input_parameters.ring import Ring
+        from blond.impedances.impedance import TotalInducedVoltage, InducedVoltageResonator
+        from blond.trackers.tracker import RingAndRFTracker
 
         gamma_transition = 1 / np.sqrt(0.00192)
         mom_compaction = 1 / gamma_transition ** 2
         circ = 6911.56
         n_rf_stations = 2
         l_per_section = circ / n_rf_stations
-        n_turns = 2
+        self.n_turns = 2
         sectionlengths = np.full(n_rf_stations, l_per_section)
         alpha_c = np.full(n_rf_stations, mom_compaction)
-        energy_program = np.array([[1*25e9,1*25e9,1.2*25e9],
-                                   [1.0*25e9,1.1*25e9,1.3*25e9]]) # todo
+        energy_program = np.array([[1 * 25e9, 1 * 25e9, 1.04 * 25e9],
+                                   [1 * 25e9, 1.01 * 25e9, 1.05 * 25e9]])
 
         ring = Ring(ring_length=sectionlengths, alpha_0=alpha_c,
-                    particle=Proton(), n_turns=n_turns,
+                    particle=Proton(), n_turns=self.n_turns,
                     n_sections=n_rf_stations, synchronous_data=energy_program)
-        self.rf_station = RFStation(ring=ring,
-                                     harmonic=[4620],
-                                     voltage=[300e6],
-                                     phi_rf_d=[0.0],
-                                     n_rf=1)
-
-        beam = Beam(ring=ring,
-                    n_macroparticles=1001,
-                    intensity=int(1e10))
-
-        self.beam = beam
-        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=2 ** 9,
+        self.rf_station = RFStation(ring=ring, harmonic=[4620], voltage=[300e6], phi_rf_d=[0.0], n_rf=1)
+
+        self.beam = Beam(ring=ring,
+                         n_macroparticles=1001,
+                         intensity=int(1e10))
+        self.n_slices = 2 ** 9
+        cut_options = CutOptions(cut_left=0, cut_right=2 * np.pi, n_slices=self.n_slices,
                                  cuts_unit='rad', rf_station=self.rf_station)
-        self.profile = Profile(beam, cut_options,
+        self.profile = Profile(self.beam, cut_options,
                                FitOptions(fit_option='gaussian'))
 
         self.R_shunt = 11897424000
-        self.f_res = 11897424000.
         self.Q_factor = 696000.0
-        self.resonator = Resonators(self.R_shunt, self.f_res, self.Q_factor)
-        print(self.Q_factor)
-        print(self.f_res)
-
-
-    def test_init(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator)
-
-    def test_init_mtw(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator,
-                                      multi_turn_wake=True)
-
-    def test_mtw_false_induced_volt(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile,
-                                      resonators=self.resonator,
-                                      multi_turn_wake=False)
-        ivr.process()
-        ivr.induced_voltage_1turn()
-        my_array = ivr.induced_voltage
-        self.assertTrue(np.any(my_array != 0.0))
-
+        fdet = -1320
+        self.f_res = 1297263703.3482404
+        self.resonator = Resonators(self.R_shunt, self.f_res + fdet, self.Q_factor)
 
     def test_total_induced_voltage(self):
-        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, resonators=self.resonator,
-                                      rf_station=self.rf_station, multi_turn_wake=False,
-                                      wake_length=self.profile.bin_size * len(self.profile.bin_centers),
-                                      array_length=2**9)
-
-        total_induced_voltage = TotalInducedVoltage(beam=self.beam, profile=self.profile,
-                                                    induced_voltage_list=[ivr], )
-        total_induced_voltage.induced_voltage_sum() # todo
-        total_induced_voltage.induced_voltage_sum()
-
-
-
-
+        from blond.trackers.tracker import RingAndRFTracker
+        buckets = 1
+        potential_min_cav = self.rf_station.phi_s[0] / self.rf_station.omega_rf[0, 0]
+        min_index = np.abs(self.profile.bin_centers - potential_min_cav).argmin()
+        self.timeArray = []
+        for turn_ind in range(self.n_turns):
+            self.timeArray = np.append(self.timeArray,
+                                       self.rf_station.t_rev[turn_ind] * turn_ind +
+                                       np.linspace(self.profile.bin_centers[0],
+                                                   self.profile.bin_centers[-1] + 2 * (
+                                                           self.profile.bin_centers[min_index] -
+                                                           self.profile.bin_centers[0]),
+                                                   self.n_slices * buckets + 2 * min_index))
+
+        ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, time_array=self.timeArray,
+                                        multi_turn_wake=True,
+                                        resonators=self.resonator,
+                                        rf_station=self.rf_station,
+                                        wake_length=len(self.timeArray) * self.profile.bin_size)
+                                        #array_length=int(len(self.timeArray) / self.n_turns))
+        ivr.induced_voltage_mtw()
 
 if __name__ == '__main__':
     unittest.main()
-- 
GitLab


From 215817b9d9ad0f38053d465e2d877831daf737a1 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 14 Mar 2025 09:26:27 +0100
Subject: [PATCH 45/48] improved test for calculation of multi-turn wake with
 single RF station only

---
 blond/impedances/impedance.py          |  1 -
 unittests/impedances/test_impedance.py | 14 +++++++-------
 2 files changed, 7 insertions(+), 8 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 937ce2e7..4fadc324 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -1076,7 +1076,6 @@ class InducedVoltageResonator(_InducedVoltage):
         self.mtw_memory[:int(self.n_time)] += self.induced_voltage
         self.induced_voltage = self.mtw_memory[:self.n_time]
 
-
     def to_gpu(self, recursive=True):
         """
         Transfer all necessary arrays to the GPU
diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index ee34e05c..4364232c 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -146,7 +146,7 @@ class TestInducedVoltageResonatorMultiTurnWake(unittest.TestCase):
         self.f_res = 1297263703
         self.harmonic = int(self.f_res * self.ring.t_rev[0])
         self.rf_station = RFStation(ring=self.ring, harmonic=[self.harmonic], voltage=[300e6], phi_rf_d=[0.0], n_rf=1)
-        self.n_slices = 2**9
+        self.n_slices = 2 ** 9
         self.beam = Beam(ring=self.ring,
                          n_macroparticles=1001,
                          intensity=int(1e10))
@@ -163,12 +163,12 @@ class TestInducedVoltageResonatorMultiTurnWake(unittest.TestCase):
         timeArray = []
         for turn_ind in range(self.n_turns):
             timeArray = np.append(timeArray,
-                                       self.rf_station.t_rev[turn_ind] * turn_ind +
-                                       np.linspace(self.profile.bin_centers[0],
-                                                   self.profile.bin_centers[-1] + 2 * (
-                                                           self.profile.bin_centers[min_index] -
-                                                           self.profile.bin_centers[0]),
-                                                   self.n_slices + 2 * min_index))
+                                  self.rf_station.t_rev[turn_ind] * turn_ind +
+                                  np.linspace(self.profile.bin_centers[0],
+                                              self.profile.bin_centers[-1] + 2 * (
+                                                      self.profile.bin_centers[min_index] -
+                                                      self.profile.bin_centers[0]),
+                                              self.n_slices + 2 * min_index))
         ivr = InducedVoltageResonator(beam=self.beam, profile=self.profile, time_array=timeArray,
                                       multi_turn_wake=True,
                                       resonators=self.resonator,
-- 
GitLab


From 88cffcc98300576da22398616cd24bd122459231 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 14 Mar 2025 09:48:19 +0100
Subject: [PATCH 46/48] changed tArray to time_array

---
 blond/impedances/impedance.py | 14 +++++++-------
 1 file changed, 7 insertions(+), 7 deletions(-)

diff --git a/blond/impedances/impedance.py b/blond/impedances/impedance.py
index 4fadc324..c328dd7f 100644
--- a/blond/impedances/impedance.py
+++ b/blond/impedances/impedance.py
@@ -978,15 +978,15 @@ class InducedVoltageResonator(_InducedVoltage):
         # If left out, the induced voltage is calculated at the times of the
         # line density.
         if time_array is None:
-            self.tArray = self.profile.bin_centers
+            self.time_array = self.profile.bin_centers
             self.atLineDensityTimes = True
         else:
-            self.tArray = time_array
+            self.time_array = time_array
             self.atLineDensityTimes = False
 
         self.array_length = array_length
 
-        self.n_time = len(self.tArray)
+        self.n_time = len(self.time_array)
 
         # Copy of the shunt impedances of the Resonators in* :math:`\Omega`
         self.R = resonators.R_S
@@ -1015,7 +1015,7 @@ class InducedVoltageResonator(_InducedVoltage):
             (self.n_time, self.profile.n_slices), dtype=bm.precision.real_t, order='C')
 
         # Call the __init__ method of the parent class [calls process()]
-        super().__init__(beam, profile,
+        super().__init__(beam=beam, profile=profile,
                          frequency_resolution=frequency_resolution,
                          wake_length=wake_length,
                          multi_turn_wake=multi_turn_wake,
@@ -1048,7 +1048,7 @@ class InducedVoltageResonator(_InducedVoltage):
                                                 self.profile.bin_centers,
                                                 self.profile.bin_size,
                                                 self.n_time, self._deltaT,
-                                                self.tArray, self._reOmegaP,
+                                                self.time_array, self._reOmegaP,
                                                 self._imOmegaP, self._Qtilde,
                                                 self.n_resonators,
                                                 self.omega_r,
@@ -1088,7 +1088,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.array(self.induced_voltage)
         self._kappa1 = cp.array(self._kappa1)
         self._deltaT = cp.array(self._deltaT)
-        self.tArray = cp.array(self.tArray)
+        self.time_array = cp.array(self.time_array)
         self._tmp_matrix = cp.array(self._tmp_matrix)
         # to make sure it will not be called again
         self._device: DeviceType = 'GPU'
@@ -1105,7 +1105,7 @@ class InducedVoltageResonator(_InducedVoltage):
         self.induced_voltage = cp.asnumpy(self.induced_voltage)
         self._kappa1 = cp.asnumpy(self._kappa1)
         self._deltaT = cp.asnumpy(self._deltaT)
-        self.tArray = cp.asnumpy(self.tArray)
+        self.time_array = cp.asnumpy(self.time_array)
         self._tmp_matrix = cp.asnumpy(self._tmp_matrix)
 
         # to make sure it will not be called again
-- 
GitLab


From 8278a96b180de617b08c4a9f5fedd98d36e816d1 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 14 Mar 2025 10:24:31 +0100
Subject: [PATCH 47/48] changed delta_E calculation to use numpy only

---
 blond/input_parameters/ring.py | 8 ++------
 1 file changed, 2 insertions(+), 6 deletions(-)

diff --git a/blond/input_parameters/ring.py b/blond/input_parameters/ring.py
index 76b3c6d7..f7137f0f 100644
--- a/blond/input_parameters/ring.py
+++ b/blond/input_parameters/ring.py
@@ -290,12 +290,8 @@ class Ring:
         else:
             # when there is more than 1 RF station, self.energy has shape (n_sections, n_turns+1)
             self.delta_E = np.zeros((n_sections, n_turns))
-            for j in range(n_sections):
-                for i in range(n_turns):
-                    if j == 0:
-                        self.delta_E[j, i] = self.energy[j, i + 1] - self.energy[-1, i]
-                    else:
-                        self.delta_E[j, i] = self.energy[j, i + 1] - self.energy[j - 1, i + 1]
+            self.delta_E[0, :] = self.energy[0, 1:n_turns + 1] - self.energy[-1, 0:n_turns]
+            self.delta_E[1:, :] = self.energy[1:, 1:n_turns + 1] - self.energy[:-1, 1:n_turns + 1]
 
         # Momentum compaction, checks, and derived slippage factors
         if ring_options.t_start is None:
-- 
GitLab


From ded95db19c6aa5865ff6f9cb7c79fcd44655b854 Mon Sep 17 00:00:00 2001
From: Elleanor Lamb <elamb@cern.ch>
Date: Fri, 14 Mar 2025 10:52:55 +0100
Subject: [PATCH 48/48] add relative file for test data

---
 unittests/impedances/test_impedance.py | 5 ++++-
 1 file changed, 4 insertions(+), 1 deletion(-)

diff --git a/unittests/impedances/test_impedance.py b/unittests/impedances/test_impedance.py
index 4364232c..e00e8c40 100644
--- a/unittests/impedances/test_impedance.py
+++ b/unittests/impedances/test_impedance.py
@@ -16,6 +16,7 @@ Unittest for impedances.impedance
 import unittest
 
 import numpy as np
+import os
 from blond.beam.profile import CutOptions, Profile, FitOptions
 from blond.impedances.impedance import InducedVoltageFreq, InducedVoltageTime, InducedVoltageResonator, \
     TotalInducedVoltage
@@ -185,7 +186,9 @@ class TestInducedVoltageResonatorMultiTurnWake(unittest.TestCase):
             longitudinal_tracker.track()
             induced_voltages.append(longitudinal_tracker.totalInducedVoltage.induced_voltage)
         induced_voltages = np.array(induced_voltages)
-        expected_voltages = np.load("./data/mtw_resonator_induced_voltage.npy")
+
+        folder = os.path.abspath(os.path.dirname(__file__)) + "/data/"
+        expected_voltages = np.load(folder + "mtw_resonator_induced_voltage.npy")
         np.testing.assert_allclose(induced_voltages, expected_voltages, rtol=1e-5, atol=1e-1)
 
 
-- 
GitLab