Model-independent particle accelerator tuning

We present a new model-independent dynamic feedback technique, rotation rate tuning, for automatically and simultaneously tuning coupled components of uncertain, complex systems. The main advantages of the method are: (1) it has the ability to handle unknown, time-varying systems, (2) it gives known bounds on parameter update rates, (3) we give an analytic proof of its convergence and its stability, and (4) it has a simple digital implementation through a control system such as the experimental physics and industrial control system (EPICS). Because this technique is model independent it may be useful as a realtime, in-hardware, feedback-based optimization scheme for uncertain and time-varying systems. In particular, it is robust enough to handle uncertainty due to coupling, thermal cycling, misalignments, and manufacturing imperfections. As a result, it may be used as a fine-tuning supplement for existing accelerator tuning/control schemes. We present multiparticle simulation results demonstrating the scheme’s ability to simultaneously adaptively adjust the set points of 22 quadrupole magnets and two rf buncher cavities in the Los Alamos Neutron Science Center (LANSCE) Linear Accelerator’s transport region, while the beam properties and rf phase shift are continuously varying. The tuning is based only on beam current readings, without knowledge of particle dynamics. We also present an outline of how to implement this general scheme in software for optimization, and in hardware for feedback-based control/ tuning, for a wide range of systems.


I. INTRODUCTION A. Motivation
It is rarely possible to build exact, deterministic input to output models for complex physical systems such as particle accelerators.It is especially difficult when the behavior of the system is influenced by many coupled parameters.Since accelerators have many coupled parameters, they are prime candidates for genetic algorithm (GA) and multiobjective genetic algorithm (MOGA) based multidimensional, nonlinear optimization schemes.In fact, MOGAs and GAs have been used to successfully optimize many aspects of particle accelerators, such as magnet and radio frequency (rf) cavity design [1], photoinjector design [2], damping ring design [3], storage ring dynamics [4], global optimization of a lattice [5], neutrino factory design [6], simultaneous optimization of beam emittance and dynamic aperture [7], and free electron laser linac drivers [8].A thorough review of GA for accelerator physics applications is given in [9].
After an accelerator design has been finalized and the accelerator has been constructed, one often encounters time-varying and nonlinear coupling effects between the imperfectly manufactured and misaligned/unknown orientation components of the accelerator.In theory, accelerator design takes a certain level of uncertainty into account.In practice however, most accelerators require postmanufacture and postinstallation tuning.This is especially the case for facilities with limited real-time diagnostics and noise measurement.In this case components may have to be retuned after each shutdown or change in operating conditions.Effects such as unknown hysteresis curves and time-varying component thermal cycling also add to system uncertainty.A particular problem faced by many accelerator systems is the arbitrary phase shift of the rf systems, a time-varying uncertainty, requiring time consuming tuning such as phase scans.The method presented here is demonstrated to automatically adapt for timevarying properties, such as phase shift.A combination of the global optimization abilities of GAs, with a local, model-independent feedback technique such as the one presented here, has the potential to improve accelerator design and performance.

B. Results of the paper
In this work, we present a simple, model-independent technique, which can aid in parameter tuning because it does not, by design, assume any particular system model for optimizing/tuning and therefore may be implemented in hardware to automatically fine-tune multiple parameters and help mitigate unmodeled disturbances and component imperfections.
For implementation, the user first defines a measurable cost function, C, to be minimized, whose analytic form may be unknown, such as the total particle loss along the length of a particle accelerator.The components p i of the vector p ¼ ðp 1 ; . . .; p m Þ are parameters by which the cost may be influenced, such as the power source current settings feeding the quadrupole magnets in the accelerator lattice.The rotation rate (RR) tuning law is Initial settings pð1Þ are chosen as usual, based on a physics model and a MOGA or other numerical optimization technique.The initial cost, C½pð1Þ, is calculated after the first run and new parameter values pð2Þ are set according to (1).
The waiting time between implementation n and n þ 1 is chosen depending on component response/settling rates, data acquisition rates, and the rate of time variation of system components due to disturbance.Intuitively, is the magnitude of a high frequency (! i ) dither being introduced into the system's dynamics, k is the gain of the controller/tuning algorithm, and Á is a time interval for digital implementation.The choices for the values of !i , Á, k, and are discussed in detail in Sec.III.We demonstrate, through a multiparticle simulation of the Los Alamos Neutron Science Center (LANSCE) low energy beam transport region, RR's ability to handle uncertainty by tuning up a 22 quadrupole lattice, and the phases of two rf buncher cavities, by minimizing a cost which is based on beam current loss along the transport region and the first two tanks of the drift tube linac, for a time-varying beam and time-varying rf phase drifts.
The RR approach may also be used as a simple numerical optimizer, in which adding new parameters p mþ1 ; p mþ2 ; . . . to (1) does not add significantly to computation time.Multiobjective optimization is implemented by replacing a single cost, C, with a combination for any number of costs, such as C 1 being total beam loss, C 2 being total transverse beam size, etc.Also, different parameters may be updated based on different rates and costs (as in Sec.IVA), where components having different sensitivities may require different values of k i and i .For example, we may implement for notational simplicity, we stick to single values of k, , and C throughout the analysis performed, which is applicable in the same manner to the above, more general scheme (2).Imposing restrictions on the parameters is straightforward and implemented as described in Sec.III.
As shown in the analysis in Sec.II B, the scheme (1) is chosen so that, with proper choices of !, k, , and Á, (1) is the finite difference approximation of which on average follows the same trajectory as the system where the convention used here and throughout the remainder of the paper is _ p ¼ @p @t .RR is related to dithering-based optimization/stabilization schemes, known in control theory as extremum seeking (ES), which have been used for optimizing unknown outputs of known, stable systems, by tuning known controllers.Originally introduced in 1922 [10], an overview of its development is available in [11].Recently, ES has been extended to perform stabilization and optimization of unknown, possible unstable systems [12].RR is a further improvement and modification of those results [13].

C. Optimization schemes
There are many existing model dependent numerical methods for multidimensional/multiparameter optimization, such as GA, MOGA, Newton-Raphson and gradient descent based on the analytic form of rC.Many optimization methods are actually built into existing accelerator design codes [14].
The main strength of RR is its model-independent nature and ability to deal with multiple parameters simultaneously, even for time-varying systems, such as thermal cycling, or unexpected component damage.Some very simple, but computationally intensive and inefficient, model-independent methods are grid and random point searches, especially for systems with many parameters.Gradient descent, based on a numerical approximation of an unknown rC is another model-independent approach, but especially in the case of a multiparameter, noisy calculation of C, may face difficulties, whereas RR is both robust to noise in C and does not need to try and estimate rC.Only samples CðnÞ are required.
Simplex fitting, in the sense that it samples many different directions in multidimensional parameter space, has the most in common with RR.A major benefit of RR is that its complexity does not grow with parameter number, regardless of the number of parameters being tuned.The scheme basically depends on three choices, the values of k, !, and .Regarding noisy data, the RR scheme is, on average, not influenced by noise, unless it happens to both match an RR parameter's perturbation frequency and be large in magnitude relative to that perturbation.Also, noise is easily handled by standard methods, such as averaging and filtering.

D. Limitations
Although RR is model independent and able to tune many parameters simultaneously, unlike GA, it is a local technique, similar to gradient descent, and may become trapped in local minimums.Therefore, we plan on exploring (in future work) a combination of GA and RR, in which a GA is first used for global optimization followed by RR for local, in-hardware tuning, to make up for modeling errors and time variation.SCHEINKER, PANG, AND RYBARCYK Phys.Rev. ST Accel.Beams 16, 102803 (2013) 102803-2

E. Organization
In Sec.II we give an explanation for the choice of the update scheme (1).In Sec.III we explain how to choose all RR parameters and describe digital implementation for general systems.In Sec.IV, we demonstrate the scheme's ability to handle uncertainty by tuning a 22 quadrupole lattice as well as two rf buncher cavities, in a simulation of the LANSCE accelerator transport region and first two tanks of the drift tube linac, for a varying input beam and drifting phase shifts.Finally, in Sec.V we provide background regarding RR analysis.

II. TUNING METHOD A. Physical motivation
It is well known that by adding a fast, small oscillation into a system's dynamics, unexpected stability properties may be achieved.The classic example is of the inverted pendulum, whose vertical equilibrium point may be stabilized by rapidly vertically oscillating the pendulum's pivot point.The dynamics of this process were first analytically described in the 1950s by Kapitza [15].The RR scheme has some similarities to this approach, in that we introduce high frequency oscillations into a system in order to force certain points of the state space to become stable equilibrium points towards which the system's trajectory converges.By abstracting this to a general state space and choosing such a point to be the minimum of a cost function, we are able to tune a wide range of systems towards various performance goals.
We start with a simple example, we do not introduce any destabilizing terms in ( 5) and ( 6), which are discussed in remark 1.To give a simple 2D overview of this method, we consider finding the minimum of a measurable function Cðx; yÞ, for which we cannot simply implement a gradient descent for the trajectory of ½xðtÞ; yðtÞ because we are unaware of its analytic form.We propose the following adaptive scheme: Note that although Cðx; yÞ enters the argument of the adaptive scheme, we do not rely on any knowledge of the analytic form of Cðx; yÞ, we simply assume that its value is available for measurement at different locations ðx; yÞ.
The velocity vector, where ðtÞ ¼ !t þ kC½xðtÞ; yðtÞ, has constant magnitude, , and therefore the trajectory ½xðtÞ; yðtÞ moves at a constant speed.However, the rate at which the direction of the trajectories' heading changes is a function of !, k, and C½xðtÞ; yðtÞ expressed as Therefore, when the trajectory is heading in the correct direction, towards a decreasing value of C½xðtÞ; yðtÞ, the term k @C @t is negative so the overall turning rate @ @t (8) is decreased.On the other hand, when the trajectory is heading in the wrong direction, towards an increasing value of C½xðtÞ; yðtÞ, the term k @C @t is positive, and the turning rate is increased.On average, the system ends up approaching the minimizing location of C½xðtÞ; yðtÞ because it spends more time moving towards it than away.
The ability of this direction-dependent turning rate scheme is apparent in the simulation of system ( 5) and ( 6), in Fig. 1.The system, starting at initial location xð0Þ ¼ 1, yð0Þ ¼ À1, is simulated for five seconds with update parameters !¼ 50, k ¼ 5, ¼ 0:5, and Cðx; yÞ ¼ x 2 þ y 2 .We compare the actual system's ( 5) and ( 6) dynamics with those of a system performing gradient descent: whose behavior our system mimics on average, with the difference max The subfigure in the bottom left shows the rotation rate, @ @t ¼ !þ @Cðx;yÞ @t , for the part of the trajectory that is bold red, which takes place during the first 0.5 seconds of simulation.The rotation of the parameters' velocity vector vðtÞ slows down when heading towards the minimum of Cðx; yÞ ¼ x 2 þ y 2 , at which time k @C @t < 0, and speeds up when heading away from the minimum, when k @C @t > 0. The system ends up spending more time heading towards and approaches the minimum of Cðx; yÞ. made arbitrarily small for any value of T, by choosing arbitrarily large values of !.The derivation of this relationship and of the rate of the gradient descent are given in Sec.V.
Towards the end of the simulation, when the system's trajectory is near the origin, Cðx; yÞ % 0, and the dynamics of ( 5) and ( 6

B. General RR scheme
For general tuning, we consider the problem of locating an extremum point of the function Cðp; tÞ: which we can measure the value of, but whose analytic form is unknown.For notational convenience, in what follows we sometimes write CðpÞ or just C instead of C½pðtÞ; t.
The explanation presented in the previous section used sinðÁÞ and cosðÁÞ functions for the x and y dynamics to give circular trajectories.The actual requirement for convergence is for an independence, in the frequency domain, of the functions used to perturb different parameters, such as sines or cosines of distinct frequencies.In what follows, replacing cosðÁÞ with sinðÁÞ throughout, or mixing sinðÁÞ and cosðÁÞ terms makes no difference.
Theorem 1.-Consider the setup shown in Fig. 2 (for maximum seeking we replace k with Àk): where The trajectory of system ( 14) approaches the minimum of Cðp; tÞ, with its trajectory arbitrarily close to that of with the distance between the two decreasing as a function of increasing !0 .Namely, for any given T 2 ½0; 1Þ, any compact set of allowable parameters p 2 K & R m , and any desired accuracy , there exists !?0 such that for all !0 > !?0 , the distance between the trajectory pðtÞ of ( 14) and " pðtÞ of ( and apply corollary 1 with respect to ! 0 and # ¼ 0:5.The trajectory of system ( 14) uniformly converges to the trajectory of where we have used the fact that mismatched terms of the form cosð! i tÞ sinð! j tÞ, 8 i, j, and terms of the form cosð! i tÞ cosð! j tÞ, and sinð! i tÞ sinð! j tÞ, 8 i Þ j weakly, uniformly converge to zero.Combining all the p i components we get Remark 1.-The stability of this scheme is verified by the fact that an addition of an unmodeled, possibly destabilizing perturbation of the form fðp; tÞ to the dynamics of _ p results in the averaged system: which may be made to approach the minimum of C, by choosing k large enough relative to the values of kðrCÞ T k and kfð " p; tÞk.Detailed stability analysis is available in [12].
Remark 2.-Although it is glossed over in the averaging analysis presented above, if one looks into the details of the proof of theorem 2, in the case of a time-varying max/min location p ? ðtÞ of Cðp; tÞ, there will be terms of the form which are made to approach zero by increasing !. Furthermore, in the analysis of the convergence of the error p e ðtÞ ¼ pðtÞ À p ? ðtÞ, there will be terms of the form Together, ( 22) and (23) imply the intuitively obvious fact that for systems whose time variation is fast, in which the minimum towards which we are descending is quickly varying, both the value of ! and of the product k must be larger than for the time-invariant case.Remark 3.-In the case of different parameters having vastly different response characteristics and sensitivities (such as when tuning both rf and magnet settings in the same scheme), the choices of k and may be specified differently for each component p i , as k i and i , without change to the above analysis.

III. GUIDELINES FOR DIGITAL IMPLEMENTATION A. Cost and constraints
The first step is to choose tunable machine parameters, p ¼ ðp 1 ; . . .; p m Þ and a cost function to be minimized, C ¼ C½p 1 ðtÞ; . . .; p m ðtÞ; t.Next, constraints for all parameters are chosen: p max ¼ ðp 1;max ;...;p m;max Þ; p min ¼ ðp 1;min ;...;p m;min Þ: Implementing initial parameter settings pð1Þ, which are chosen based on the physics model and numerical methods, allows one to calculate C½pð1Þ.The iterative update scheme is then which is based on the finite difference approximation of the derivative: which, according theorem 1 will drive the system towards a minimum of C. The constraints are implemented by checking the updated parameters at each step and confining them to their bounds if necessary: It is important that !i ) kC, so that the adaptive scheme is operating on a faster time scale then and able to adapt to time variation of the cost function.Because RR depends on distinguishing between different frequency components of the cost, the !i should be chosen within such bounds such that it is possible to implement a small enough Á obeying ensuring that at least 20 iterations (10 Â the Nyquist sampling rate) are required to perform one complete cosine oscillationintheiterativescheme(24).Choosingsmallervalues of Á results in smoother parameter oscillation and more iterativestepsrequiredforconvergence,largervaluesofÁ speedup the convergence, but may destabilize the overall scheme.
According to theorem 1, the only requirement on the choices of !i is that they are big enough and distinct, but in practice, the more harmonically independent they are (such as !i Þ 2! j for all i Þ j) the better.The sensitivity to frequency independence is different for every system and depends on the coupling between different components.One simple method is to choose a scaling factor, !0 , and !i ¼ !0 r i , where the values r i are distinct.
The update scheme (24) is only valid as a finite difference approximation as in (25) if Á ( 1 and !i ) ffiffiffi ffi p ; kC.Taking into account (26), we choose a large value of !0 , relative to kC, typically at least is a safe choice, where C½pð1Þ is the initial cost calculated based on initial parameter settings pð1Þ.
The rate of convergence is proportional to the product k, increasing either k or speeds up convergence, as long as they are not too big relative to the value of !0 , so that the finite difference is an accurate approximation of the derivative.If, after !0 has been chosen, the convergence is too slow, or if a local minimum is suspected, k or may be increased, with the possible need to increase !0 as well.The vector p is moving through the parameter space R m in ellipses with approximate major axes of magnitude ffiffiffi !p , increasing causes larger steady state parameter oscillations, which is not a problem if the adaptation is turned off following successful convergence.
The choices described above may vary from system to system based on sensitivity and initially may be an iterative process.A good approach is to first fix values of k and , define the various relationships (26) and ( 27), and increase !0 if necessary until the scheme is stable.Once convergence begins, if it is too slow, or if the cost is not sensitive enough, k and may be increased, with a possible necessity to increase !0 as well to maintain stability.

C. Digital resolution
Although the analytic form of CðnÞ may be unknown, at each iteration the parameters are perturbed by quantities with known bounds: 2 n b , then, as cosð Þ varies between 0 and 1, it is possible for the parameter value to take N discrete steps of minimum resolution M i 2 n b .

D. Normalization of parameters
Different parameters p i may require individual values of k i and i , in which case normalizing the parameters to within ½À1; 1 bounds may be useful.For example, at each step n, one may compute the cost CðnÞ based on parameter settings pðnÞ, then translate into the scaled parameters p s ðnÞ: where C p;i ¼ p i;max þp i;min 2 and D p;i ¼ p i;max À p i;min , bounding each parameter within ½À1; 1.We then perform the RR update, force the scaled parameters to satisfy the constraints À1 and 1, and transform back into unscaled parameter values in order to calculate the cost for the next iteration:

IV. TUNING 22 QUADRUPOLE MAGNETS AND 2 BUNCHER CAVITIES
In this section we present simulation results of using the RR scheme to tune up the 22 quadrupole magnets and two buncher cavities in the Los Alamos linear accelerator H þ transport region, a simplified schematic of which is shown in Fig. 3.The simulations were done using a graphics processing unit-accelerated online beam dynamics simulator [16,17], which is being developed to predict beam properties along the linac using real-time machine parameters.It can serve as a virtual beam experiment environment and contribute to the cost being minimized by the RR optimizer, by providing pseudo-real-time estimates of beam sizes and current information in parts of the machine where diagnostics are not available.Currently being demonstrated on the LANSCE low energy beam transport (LEBT) and drift tube linac (DTL), simulating a bunch of 32K macroparticles through the LEBT or DTL takes fractions of a second, which is 40 times faster than the simple CPU version of the code.

A. Magnet tuning for beam transport
In a first, simple demonstration of the technique, we perform a simulation of only the LEBT, with all initial magnet current set points set to 0A, and allowed to tune up based purely on the RR scheme as described above, in which the four costs (j ¼ 1; 2; 3; 4) being minimized, C j ¼ ðI j À 0:013Þ 2 ; (32) were the square of the difference between initial beam current 0:013A and total current making it through various parts of the transport region, at which diagnostics are available.With reference to Fig. 3, the current is sampled at four locations, I 1 , following Q 6 , I 2 following Q 10 , I 3 following Q 18 , and I 4 at the end of the transport region.The magnets (i ¼ 1; . . .; 22) were then updated according to where and S i ¼ C 4 for Q 19 -Q 22 , so that magnets only saw costs which they were able to influence.For the tuning parameters, we chose k ¼ 250000, so that the amplified costs kS j in (33) took values between 0 and 300.The ! i were chosen kC max > 20. Figure 4 shows the evolution of the surviving beam current at the end of the transport region during the RR tuning scheme.Figure 5 shows the evolution of the magnet current inputs.Figure 6 shows the rms beam size through various parts of the transport region at the end of RR tuning, and Fig. 7 compares the RR found magnet settings to that of the tune up in 2011.
This example demonstrates some of the strengths and limitations of the scheme, and the importance of cost function choice.Although the cost has been minimized and almost all current is making it to the end of the transport region, the beam is beginning to diverge and in this form would not be matched to the DTL following the transport region.In practice it is of course better to start with physics-model based initial parameters; this simulation was conducted starting with all magnet settings at zero in order to fairly demonstrate the model-independent abilities of the RR scheme.The next simulations start with the 2011 tune up for the magnet settings and use current monitors following two tanks of the DTL, in which case surviving beam corresponds with well-matched beam.

B. Magnet and rf buncher cavity tuning
To demonstrate the use of this scheme for fine-tuning of machine settings, we used machine settings found during the 2011 tune up procedure, but with a slightly different beam and incorrectly phased buncher cavities.The magnets were initialized to the values recorded from one of the 2011 machine turn on tuning periods.We set the phase settings for the buncher and prebuncher to zero, which typically must be retuned at each turn on, by a phase scan, to take care of arbitrary phase shift.
We used only the surviving current at the end of the second tank of the drift tube linac to create our cost; our tuning procedure for the parameters was where i ¼ m for the magnets and i ¼ b for the buncher phases.In both cases CðnÞ ¼ ðI end À 15 mAÞ 2 : For the tuning parameters, we chose k ¼ 605000, m ¼ 25, b ¼ 550.The ! i were chosen as !0 r i , with !0 ¼ 2000 and r i uniformly distributed between 2.5 and 4.3, Á ¼ 2% 20! 24 .With these values, !min kC max > 35.With an initial beam current of 15 mA, the typical surviving current after machine tune up is roughly 80% or 12 mA.After 2000 simultaneous iterations on these 24 parameters (22 quads, 2 buncher phases), the surviving current at the end of tank 2 was 12.25 mA.The results of the optimization procedure are shown in Figs.8-12.From Figs. 9 and 10 we see that only minor adjustments are made to magnet settings compared to the rf phases.Figure 11 shows that the transverse beam size has further focused throughout the transport region and the transverse match to the DTL has slightly improved.Figure 12 compares surviving beam current at the end of tank 2 of the DTL before and after tuning.
Step Number Step Number  In order to demonstrate the adaptive tuning abilities of the scheme, we started with matched beam settings and varied both the characteristics of the input beam and added a time-varying phase drift to each buncher cavity.
Figure 13 shows the initial and final beam properties at the entrance to the transport region, during which RR adaptive tuning maintains beam focus and matching.Figure 14 shows the phase shift of the bunchers with and without tuning.These changes took place starting at step 1000 and finished at step 19000, with beam properties staying constant before and after the interval.Also, during this beam changing process, the phase of the first buncher was made to drift by 30 deg and that of the second by 35 deg, as seen in Fig. 14.

Initial Machine Settings
After Optimization 80% of Initial Current   The drift of beam characteristics and buncher phase shifts took place over 18000 time steps, which for a conservative magnet/phase update rate of 1 Hz translates into drastically changing accelerator and beam properties over the course of just 5 hours.All tuning parameters were maintained exactly the same as in the previous example.
Figure 15 shows the evolution of the magnet gradients throughout the process and Fig. 16 compares the initial and final beam profiles.In Fig. 17 we see that adaptive RR tuning is able to maintain $12 mA of surviving beam during the time-varying beam and phase, whereas almost all of the beam is lost without tuning.

V. ANALYTIC BACKGROUND
We briefly recall the functional analysis result of Kurzweil and Jarnik [18], which allows one to relate the trajectories of a highly oscillatory system to those of a simplified Lie bracket averaged system.

VI. CONCLUSIONS AND FUTURE WORK
Because of the global optimization ability of MOGAs and the simple, fast, abilities of RR, we think that a combination of RR and MOGA techniques can be a very powerful numerical optimization method, in which once the neighborhood of a global optimal solution has been determined by MOGA, RR may be implemented for fine-tuning, and finally in hardware to compensate for unmodeled system characteristics.In the future, following upgrades to the LANSCE digital control system and network, we plan to test this algorithm on various accelerator components.
FIG.1.The subfigure in the bottom left shows the rotation rate, @ @t ¼ !þ @Cðx;yÞ @t , for the part of the trajectory that is bold red, which takes place during the first 0.5 seconds of simulation.The rotation of the parameters' velocity vector vðtÞ slows down when heading towards the minimum of Cðx; yÞ ¼ x 2 þ y 2 , at which time k @C @t < 0, and speeds up when heading away from the minimum, when k @C @t > 0. The system ends up spending more time heading towards and approaches the minimum of Cðx; yÞ.

FIG. 4 .
FIG.4.The surviving current at the end of the beam transport over 2500 iteration steps is shown for an initial beam current of 13 mA.

FIG. 5 . 2 FIG. 6 .
FIG.5.Evolution of the magnet current settings to the magnets over 2500 iteration steps.

FIG. 7 . 7 C
FIG.7.Magnet settings at the end of the iterative tuning scheme compared to 2011 tune up settings.
FIG.15.Evolution of the magnet gradients before, during, and after variation of beam and phase delay parameters.

12 FIG. 16 .
FIG.16.Comparison of rms beam size along the accelerator for the 2011 tune-based magnet settings (dashed) and after the beam initial conditions have changed and RR tuning has refocused and matched the beam (solid).
of or large values of !.A detailed overview of how to choose the values k, , and ! is given in Sec.III.Convergence towards a maximum, rather than a minimum, is achieved by replacing k with Àk.
FIG. 2. Tuning of the ith componentp i of p ¼ ðp 1 ; . ..; p n Þ 2 R n .The symbol 1s denotes the Laplace transform of an integrator, so that in the above diagram p i ðtÞ ¼ p i ð0Þ þ R t 0 u i ð(Þd(.