# Low EffortLi Nuclear Fusion Plasma Control Using Model Predictive Control Laws

Authors: Izaskun Garrido; Aitor J. Garrido; Jesús A. Romero; Edorta Carrascal; Goretti Sevillano-Berasategui; Oscar Barambones
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

One of the main problems of fusion energy is to achieve longer pulse duration by avoiding the premature reaction decay due to plasma instabilities. The control of the plasma inductance arises as an essential tool for the successful operation of tokamak fusion reactors in order to overcome stability issues as well as the new challenges specific to advanced scenarios operation. In this sense, given that advanced tokamaks will suffer from limited power available from noninductive current drive actuators, the transformer primary coil could assist in reducing the power requirements of the noninductive current drive sources needed for current profile control. Therefore, tokamak operation may benefit from advanced control laws beyond the traditionally used PID schemes by reducing instabilities while guaranteeing the tokamak integrity. In this paper, a novel model predictive control (MPC) scheme has been developed and successfully employed to optimize both current and internal inductance of the plasma, which influences the L-H transition timing, the density peaking, and pedestal pressure. Results show that the internal inductance and current profiles can be adequately controlled while maintaining the minimal control action required in tokamak operation.


---

## Body

## 1. Introduction

The economy of future tokamaks [1] may have to rely on having large bootstrap current fractions and/or pulsed operation [2], with limited power available from noninductive current drive actuators [3]. However, the transformer primary coil can assist in reducing the power requirements of the noninductive current drive sources needed for current profile control, since the general shape of the current profile can be easily transiently manipulated by the transformer [4, 5].Thus, a natural extension of the existing magnetic control systems (plasma current, shape and position control using the PF coil system [6–10]) is to transiently add the control of the magnetic field structure inside the plasma using the transformer primary coil. The scheme will be particularly useful for the startup and termination phases of future pulsed reactors such as ITER [11, 12], or ignition designs relying on a fast plasma current ramp, such as Ignitor, as well as of present day tokamak research facilities.To test these ideas, the system should use the flux change provided by the transformer primary coil to control either the electric current or the internal inductance of the transformer’s secondary plasma circuit load. The internal inductance control [13] is used to regulate the slow flux penetration due to the skin effect, providing first-order control over the shape of the plasma current density profile in the highly conductive plasma. In this context, this article implements model predictive control techniques using a state space model, with the goal to create the necessary control instruments for ITER scenario operation that complements the great effort that is being deployed within advanced scenario modeling.In particular, new control models have been developed in this research, which will help to improve the behaviour of the major characteristics with subsequent plasma stabilization in agreement with the objective of this research, which is to regulate throughout the pulse the internal inductance and the plasma current sensitivity of the advanced scenarios with respect to transport models and physical assumption. Having into account that the promising results that have been obtained for the proposed state space system based on lumped parameter simulations; this encouraging outcome reflects the importance of internal inductance control using the OH coil for ITER like advanced scenarios, due to the limited power available for noninductive current drive actuators. Thus, more complexity could be added to the space state system in order to accurately model the relation between internal inductance and current profiles and, to a lesser extent relative to the project aim, sensitivity to the L-H transition timing, to the density peaking and pedestal pressure. This research is even more relevant if one considers the project timing in the current ITER development timeline when the need for control emerges as a key factor in the development of the Missions established in the roadmap to fusion.(i)
Mission 1: Plasma Regimes of Operation. Plasma regimes of operation-based on the tokamak configuration- for reactor application need to achieve high fusion gain by controlling plasma instabilities and minimizing the energy losses due to small-scale turbulence [14–23].(ii)
Mission 2: A Reliable Solution to the Problem of Heat Exhaust in the Fusion Power Plant. Achieving conditions, in which a net surplus of fusion energy is produced, requires control techniques so as to maintain plasmas at high density and temperature for a few hours or even in steady state.(iii)
Mission 3: DEMO Design. Specific activities to demonstrate the control of plasma regimes of operations with DEMO relevant systems, such as plasma control using the transformer primary in order to reduce the number of actuators.Control of the plasma inductance is an essential tool for current profile control in any tokamak reactor, which can be used to extend the pulse duration, access to advanced regimes, reduce the growth rate of vertical instability, and improve the reproducibility of the experiment. In particular, it is paramount in ITER for the following reasons:(i)
The vertical stability should be provided largely by passive means, mainly through induction of image currents that occurs inserting small quantities of material with low resistivity in the steel structure of the vacuum chamber. In addition, it must be actively controlled, mainly by external coils that generate a radial magnetic field component but also by faster internal coils unfiltered by the vacuum vessel. The growth rate of the vertical instability is a positive function of the internal inductance. In order to increase the stability margins, it is therefore essential to keep the inductance as low as possible in the transitional phases, especially during the plasma current ramp-down phase of the discharge.(ii)
The lower the inductance of a given electrical circuit is, the smaller the transient voltage that must be applied to reach a given current in a given time is. To save volt seconds in the tokamak transformer primary coil, it is therefore necessary to reduce the plasma inductance. Due to the current diffusion processes, there exists a correlation between resistive and inductive flux consumption, which are related through the Ejima coefficient, so that lowering the resistive flux (for instance by plasma heating during the ramp up) also lowers the inductive flux consumption. Hence, the internal inductance of the resulting plasma is reduced and more flux is available to extend the flat-top phase of the discharge.(iii)
Broad current profiles are required to access improved confinement regimes withq>1 in the central region of the plasma, which correspond with low internal inductance plasmas. However, access to these regimes requires large flux consumption during the ramp up phase, shortening the available flux to execute the flat-top phase of the discharge. Thus, in order to achieve low inductance the existing plants require a combination of plasma heating during the ramp-up and fine-tuning of the profile with noninductive current drive in the flat-top phase.The importance of the inductance control is backed-up by the quantity and quality of scientific research performed in the last years [12, 24–27]. The existing correlation between plasma current ramp rates and inductance changes suggest the idea of using the plasma current ramp rate as a virtual actuator that serves to bring the system to the desired regime of operation. It is expected that the use of the CS or OH coil significantly reduces the power requirements of the noninductive current drive sources required for current profile control, since the low order moments of the current profile can be easily manipulated by the transformer, which in turn, avoids the decrement in plasma confinement when increasing total injected power [28].


## 2. System Description

A tokamak is a fusion reactor with a toroidal chamber where the current is induced by coils acting as the primary circuit of a transformer while the plasma itself is the secondary circuit. The magnetic field that confines this plasma is created in toroidal direction by coils located along the torus (toroidal field coil) together with another field perpendicular to the first one, created mainly by the plasma current (poloidal field). Thus, the resulting magnetic field lines are composed by the combination of these two fields (poloidal and toroidal) and present a helical shape along the torus, so that the particles pass alternately by internal and external areas of the torus.This section includes a brief description of the numerical model that will be used to control the plasma internal inductance. The model considers the lump parameter formulation described in [27], that is, from an electrical point of view the tokamak is modeled as a toroidal transformer primary coupled with the plasma ring by a mutual inductance M and the plasma acts as the secondary RL circuit, where R and L denote the plasma resistance and inductance, respectively [29–32]. The plasma is maintained using poloidal field discharges with the particularity that the total inductance consists of the sum of a constant term (external inductance) due to the inductance of transmission lines and wiring which are geometrical factors, plus a variable term corresponding to the internal inductance of the system.Besides, the shape and position of the current at a given time determine the value of the internal inductance, which is a measure of the width of the current profile. Therefore, setting correctly the internal inductance values adjusts the kinematics and shape of the current and the potential drops in the electrical circuit. This task may be improved by using modern advanced control laws already used in nuclear fusion and other energy related areas -mainly SMC robust and MPC controllers-, beyond the traditionally employed PID schemes [33–41]. Since the internal inductance has a direct bearing on the stability at a given equilibrium, there is a growing interest to exploit internal inductance control as a means to extend the duration of tokamak plasma discharges [42], to reduce the growth rate of the vertical instability of elongated plasmas [12, 43, 44] and to guarantee access to advanced tokamak scenarios [45].The first step when developing a controller is to find a suitable mathematical model of the system. The numerical model under study is the lumped parameter model for the internal inductance of plasma current that was published in [29], which is derived considering energy conservation and flux balance together with a first order approximation for the dynamics of the flux diffusion and has been validated with experimental data from TCV tokamak. Let’s recall the equations used to model the plasma current and its internal inductance evolution, and introduce the state space vector x=(x1,x2,x3,x4)T with
(1)x1=Li,x2=Ip,x3=Vc,x4=∂Vc∂t,
where Li denotes the internal inductance, Ip the total plasma current, Vc is the variation of the equilibrium flux Vc=-∂ψc/∂t and ψc is the weighted flux average for the current density enclosed by the plasma boundary Ω(2)ψc=∫ψidSI.
The following lumped state space system model is derived
(3)x˙=fx+gxuVB,R,I^,y=LiIpT
given that VB is the boundary loop voltage, R is the plasma resistance, and I^ the noninductive current
(4)fx=2Rx2-I^-x4x2x4-2Rx2-I^Le+x1-Rx2-I^-x4ω2k-1x3+kx1x2-ω2Rx2-I^α-βx4,(5)α=TRk-1-2LeTBkLe+x1,β=2δω+2LeTBkω2Le+x1,(6)gx=0-1Le+x102LeTBkω2Le+x1T
denoting Le the external inductance. This system is accurate, except for an ad hoc first order approximation for the flux diffusion dynamics by considering the voltage at the equilibrium surface as a function of τ and k, which represent, respectively, the time constant and the gain of the system
(7)Vc=x3-kx1x2τ+VB.


## 3. MPC Scheme

The general design objective of the Model Predictive Control is to compute a trajectory of a future manipulated variable, the boundary loop voltage,VB, to optimize the future behavior of the internal inductance. The optimization is performed within a limited time window by giving the initial conditions at the beginning of the time window. Although the optimal trajectory of future control signal, VB, is completely described within the moving horizon window, the actual control input to the plasma only takes the first sample of the plasma current, while neglecting the rest of the trajectory.In order to make the best decision, a criterion is needed to reflect the objective. The objective is related to an error function based on the difference between the desired and actual responses. This objective function is often called the cost function,J, and the optimal control action is found by minimizing this cost function within the optimization window.The model predictive control system is designed based on the mathematical model of the tokamak, (3). The model to be used in the control system design is a discretization of the linearized system. By using a state-space model, the current information required for predicting ahead is represented by the state variable at the current time.Considering a zero order hold discretization, the system (4) may be rewritten in matrix form as
(8)xdk+1=Adxdk+Bdudk,yk=Cdxdk,
where u is the manipulated variable that correspond to the plasma current, y is the internal inductance and xd is the space state vector. Since this model has ud(k) as input, it will be changed to suit a design purpose in which an integrator is embedded. The difference of the state space equation is
(9)Δxdk+1=AdΔxdk+BdΔudk,
where Δxd(k) denotes the difference of the state space vector Δxd(k+1)=xd(k+1)-xd(k) and Δud(k) the difference of the control variable Δud(k)=ud(k)-ud(k-1).Therefore, the input to the state space model isΔud(k). To connect the increment in the state variables to the output y(k), a new state space vector is xk=Δxd(k)Ty(k)Tand the output vector is
(10)yk+1-yk=Cdxdk+1-xdk=CmΔxmk+1=CdAdΔxdk+CdBdΔudk.
Thus, system (8) may be rewritten as the so called augmented system model
(11)Δxdk+1yk+1=Ad0CdAd1Δxdkyk+BdCdBdΔudk,yk=01Δxdkyk.
Given this mathematical model, the general expression takes the form of the incremental differential matrix equation
(12)xk+1=Axk+Buk,yk=Cxk.
The next step is to calculate the predicted internal inductance output with the future boundary loop voltage, VB, control signal as the adjustable variable. This prediction is described within an optimization window, which may be denoted Np.


### 3.1. Prediction of the Internal Inductance

Given the initial values for the state space vector at a given sampling instant, it is necessary to determine the future control trajectory(13)uk,uk+1,…,uk+Nc-1,
where Nc, the control horizon, determines the number of parameters used to predict the state vectors within the given optimization window
(14)x(k+1∣k),x(k+2∣k),…,x(k+Np∣k)
being the control horizon Nc smaller than the prediction horizon Np.Iterative substitution on the state space model (12) leads to
(15)x(k+Np∣k)=ANpxk+ANp-1Buk+ANp-2Buk+1+⋯+ANp-NcBuk+Nc-1
and the corresponding system output is
(16)yk+Npk=Cxk+Npk.
Besides, the predicted internal inductance is formulated in terms of the space vector initial conditions and the plasma current at the control horizon. Therefore, (12) may be rewritten in compact form as
(17)y^=Fxk+Gu^k,
where
(18)F=CACA2CA3⋯CANpT;G=CB00⋯0CABCB0⋯0CA2BCABCB⋯0⋮CANp-1BCANp-2BCANp-2B⋯CANp-NcB.


### 3.2. Optimization

For a given referencerk at a given time, within a prediction horizon the control system aims to bring the internal inductance as close as possible to the reference signal
(19)r^=rk⋯rk+Np-1T
and the corresponding cost function
(20)J=r^-y^Tr^-y^+u^TRu^
defines the objective of minimizing the output error and the control action. However, this minimization might be weighted by the symmetric positive definite matrix R. When R has entries with small value, it means that there is no penalty on the plasma current control action but the objective is to keep the internal inductance error as low as possible all along the prediction horizon.The cost function (20) that will determine the optimal control action in (17) may be computed as the u` that will minimize the cost function J given by
(21)J=r^-FxTr^-Fx-2u^TGTr^-Fx+u^TGTG+Ru^,
where the minimum is obtained when the first derivative with respect to the control
(22)∂J∂u`=-2GTr^-Fx+2GTG+Ru^
equals zero, which happens when
(23)u^=GTG+R-1GTr^-Fx
being (GTG+R)-1 the hessian matrix in the optimization.


## 3.1. Prediction of the Internal Inductance

Given the initial values for the state space vector at a given sampling instant, it is necessary to determine the future control trajectory(13)uk,uk+1,…,uk+Nc-1,
where Nc, the control horizon, determines the number of parameters used to predict the state vectors within the given optimization window
(14)x(k+1∣k),x(k+2∣k),…,x(k+Np∣k)
being the control horizon Nc smaller than the prediction horizon Np.Iterative substitution on the state space model (12) leads to
(15)x(k+Np∣k)=ANpxk+ANp-1Buk+ANp-2Buk+1+⋯+ANp-NcBuk+Nc-1
and the corresponding system output is
(16)yk+Npk=Cxk+Npk.
Besides, the predicted internal inductance is formulated in terms of the space vector initial conditions and the plasma current at the control horizon. Therefore, (12) may be rewritten in compact form as
(17)y^=Fxk+Gu^k,
where
(18)F=CACA2CA3⋯CANpT;G=CB00⋯0CABCB0⋯0CA2BCABCB⋯0⋮CANp-1BCANp-2BCANp-2B⋯CANp-NcB.


## 3.2. Optimization

For a given referencerk at a given time, within a prediction horizon the control system aims to bring the internal inductance as close as possible to the reference signal
(19)r^=rk⋯rk+Np-1T
and the corresponding cost function
(20)J=r^-y^Tr^-y^+u^TRu^
defines the objective of minimizing the output error and the control action. However, this minimization might be weighted by the symmetric positive definite matrix R. When R has entries with small value, it means that there is no penalty on the plasma current control action but the objective is to keep the internal inductance error as low as possible all along the prediction horizon.The cost function (20) that will determine the optimal control action in (17) may be computed as the u` that will minimize the cost function J given by
(21)J=r^-FxTr^-Fx-2u^TGTr^-Fx+u^TGTG+Ru^,
where the minimum is obtained when the first derivative with respect to the control
(22)∂J∂u`=-2GTr^-Fx+2GTG+Ru^
equals zero, which happens when
(23)u^=GTG+R-1GTr^-Fx
being (GTG+R)-1 the hessian matrix in the optimization.


## 4. Results

Considering that tokamaks are nonlinear systems, with limited power available from noninductive current drive actuators, their operation would benefit from plasma with internal inductance as small as possible. Thus, it was decided to test the lump parameter model with a predictive control scheme able to produce a controlled ramp up with optimal internal inductance value. In particular, predictive mode control has been considered as an efficient method for dealing with plasma control problems in tokamaks, due to its ability to ensure that both the control action and the solution error are minimal.In order to make the best decision, a criterion is needed to reflect the relations between the control effort and an error function based on the difference between the desired and actual responses. This objective or cost function,J, determines the optimal control action within the optimization window. Besides, increasing the control horizon will smooth the control strategy, taking it longer to reach the stationary.Therefore, the optimization window and entries in the penalty matrixR will be considered as the main criteria relevant to the control actions. The given reference rk is assumed to be a step that jumps from 0.7 0,7 μH to 1,4 μH at k=25.


### 4.1. Increasing the Size of the Optimization Window

The simulation results represented in Figure1 show the time evolution for the desired internal inductance that is obtained using the proposed MPC. It can be appreciated that after a transitory time the internal inductance tracks the desired reference for a relatively small prediction horizon, which means that the computational effort is kept low. Even more, it may be observed by comparison between Figures 1 and 2 that increasing the prediction horizon might also bring unwanted side effects like overshooting and higher control effort.Internal inductance (a), control (b)Np=5.
(a)
(b)Internal inductance (a), control (b)Np=30.
(a)
(b)


### 4.2. Increasing the Penalty in the Control Action

In this subsection we perform a study to test the advantages of penalizing the control action. It may be seen in Figure3 that a strong penalization is not satisfactory because the internal inductance is not properly controlled. However, comparing Figures 4 and 5 it is seen that some penalization is beneficial because the closed loop preserves the control on the internal inductance while saving much needed power in the control action.Internal inductance (a), control (b) penalty 1.
(a)
(b)Internal inductance (a), control (b) penalty 0.1.
(a)
(b)Internal inductance (a), control (b) penalty 0.001.
(a)
(b)


### 4.3. Comparison with Traditional PID-Based Controllers

In this subsection simulation results using a modified traditional PID-based controller with approximate derivative action are presented for comparison purposes. The method used for the PID tuning is explained in [46]. In particular, the Ziegler-Nichols frequency response method has been employed. Using this procedure, some initial values for the controller parameters were obtained and later on refined in an experimental trial and error basis. This well-known method offers a distinguished advantage since it can be performed from an external point of view, eliminating the use of an exhaustive internal modeling. In this sense, it is worth to revise the work of Professor Åström, the author of the best-selling trilogy on PID control, where, after an exhaustive analysis of the existing methodologies, it was concluded that “there are already many tuning methods available, but a replacement of the Ziegler-Nichols method is long overdue.”In this context, Figure6 presents the same internal inductance control case implemented in the previous sections, where the controller has been tuned as indicated above to solve the trajectory tracking problem in order to eliminate the steady state error. However, it may be observed from Figure 6 that the control action is much more expensive than in the cases of the MCP previously implemented. Besides, the trajectory tracking worsens when a lower control action is considered. In contrast, MPC intrinsically ensures an adequate trajectory tracking while minimizing the control signal, providing an optimal solution.Internal inductance (a), control (b) PID-based.
(a)
(b)


## 4.1. Increasing the Size of the Optimization Window

The simulation results represented in Figure1 show the time evolution for the desired internal inductance that is obtained using the proposed MPC. It can be appreciated that after a transitory time the internal inductance tracks the desired reference for a relatively small prediction horizon, which means that the computational effort is kept low. Even more, it may be observed by comparison between Figures 1 and 2 that increasing the prediction horizon might also bring unwanted side effects like overshooting and higher control effort.Internal inductance (a), control (b)Np=5.
(a)
(b)Internal inductance (a), control (b)Np=30.
(a)
(b)


## 4.2. Increasing the Penalty in the Control Action

In this subsection we perform a study to test the advantages of penalizing the control action. It may be seen in Figure3 that a strong penalization is not satisfactory because the internal inductance is not properly controlled. However, comparing Figures 4 and 5 it is seen that some penalization is beneficial because the closed loop preserves the control on the internal inductance while saving much needed power in the control action.Internal inductance (a), control (b) penalty 1.
(a)
(b)Internal inductance (a), control (b) penalty 0.1.
(a)
(b)Internal inductance (a), control (b) penalty 0.001.
(a)
(b)


## 4.3. Comparison with Traditional PID-Based Controllers

In this subsection simulation results using a modified traditional PID-based controller with approximate derivative action are presented for comparison purposes. The method used for the PID tuning is explained in [46]. In particular, the Ziegler-Nichols frequency response method has been employed. Using this procedure, some initial values for the controller parameters were obtained and later on refined in an experimental trial and error basis. This well-known method offers a distinguished advantage since it can be performed from an external point of view, eliminating the use of an exhaustive internal modeling. In this sense, it is worth to revise the work of Professor Åström, the author of the best-selling trilogy on PID control, where, after an exhaustive analysis of the existing methodologies, it was concluded that “there are already many tuning methods available, but a replacement of the Ziegler-Nichols method is long overdue.”In this context, Figure6 presents the same internal inductance control case implemented in the previous sections, where the controller has been tuned as indicated above to solve the trajectory tracking problem in order to eliminate the steady state error. However, it may be observed from Figure 6 that the control action is much more expensive than in the cases of the MCP previously implemented. Besides, the trajectory tracking worsens when a lower control action is considered. In contrast, MPC intrinsically ensures an adequate trajectory tracking while minimizing the control signal, providing an optimal solution.Internal inductance (a), control (b) PID-based.
(a)
(b)


## 5. Conclusions

The need for an optimal and robust control emerges as a key factor in the development of future fusion reactors as ITER. In this paper, a tokamak transformer model including a lumped parameter formulation for the skin effect has been used for the design of a novel MPC scheme that uses the OH coil current ramp rate to control the current and internal inductance of the plasma. Internal inductance is used to regulate the slow flux penetration in the secondary plasma due to the skin effect process, providing some basic form of current profile control. The excellent results obtained with this MPC approach are very encouraging since MPC naturally fits the need for an optimal internal inductance coupled with minimal control action. It is shown that MPC algorithms are superior compared to other classical control methods due to its ability to ensure that both the control action and the solution error are minimal which is particularly useful in the advanced tokamaks case, where limited power is available. In this sense, it is also shown that in the case of the MPC-scheme it is possible to use a relatively small control horizon while preserving the optimal behavior.



---

*Source: 527420-2015-02-23.xml*