A practical protocol to emulate a reactor scenario on present machines, with application to the ASDEX Upgrade tokamak via predictive modeling

In this work, a novel practical strategy to emulate a reactor scenario on present tokamak experiments is presented. A recipe how to scale several relevant parameters from a hypothetical reactor scenario down to present devices is discussed. Equivalence between the energy flux channels is detailed, and the practical actuation scheme is presented. The application of the proposed protocol on the ASDEX Upgrade tokamak is shown foremost using the virtual flight simulator Fenix, with practical experiments planned for future campaigns.


Introduction
Recently there has been much advance on our capabilities to emulate yet unobservable systems employing several types of similarity/duality arguments. One of the most outstanding examples is the successful experiment on the a See Stroth et al 2022 (https://doi.org/10.1088/1741-4326/ac207f) for the ASDEX Upgrade Team. * Author to whom any correspondence should be addressed.
Original Content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
Sachdev-Ye-Kitaev (SYK) model-transversable wormhole duality protocol that has been performed on quantum computers at Google [1,2], where a set of q-bits, divided into two entangled sub-sets on the two sides of the system, is able to transfer (teleport) information from a q-bit coupled from one side, to another receiving q-bit on the other side. In the field of nuclear fusion research, similarity, or identity, arguments [3,4] paved the way to systematically compare the confinement and performance of different experiments. These similarity arguments have also been employed to perform studies on isotope physics [5], momentum transport [6], for example. While in the quantum experiment case, the issue was to emulate a completely different physical system which happens to be represented by the same set of equations, in the fusion case one wants to emulate the same physical system but on a smaller scale (eventually having to live with compromises on different normalized parameters as well).
Nevertheless, both worlds are accomunated by the fact that scalings and limits, approximations and valid physical arguments are used to bring the two systems on an overlapping trajectory of 'describing equations/relevant normalized parameters'. This process of reduction of the initial problem, which cannot be reproduced in present experiments, to something which instead can be executed on our devices, involves several first principles and empirical arguments. In this work, we address these specifically for the goal of replicating a reactor scenario on a present existing device, as much as possible.
The structure of the present work is as follows: in section 2, the general 'reactor emulation protocol' (REP) is described. In section 3, application to the ASDEX Upgrade tokamak (AUG) parameters is discussed, upon using the flight simulator Fenix to carry out the virtual experiments. Finally, in section 4, conclusions are drawn.

The reactor emulation protocol (REP)
It is known that a physically sound way to obtain two equivalent plasmas (in terms of the physics) on different devices, or even in the same device but with different engineering parameters, is to employ similarity arguments [3]. We first review here briefly the classical argument in a slightly different form, which is useful for the following development. Notice that in the following we assume the aspect ratio A = R/a (i.e. ratio of major to minor plasma radius) to be a constant number (typically A ≈ 3 for standard tokamak devices). As such, A will not appear in the formulas shown (but it is straightforward to include it). Moreover, a fusion reactor will operate in DT mixture fuel, whereas we run the experiments in D. Isotoperelated effects and appearance of main ion mass in the formulas will be ignored. However, similar to A, isotope dependencies can be re-introduced in a straightforward manner at least qualitatively.

The standard similarity argument
As said previously, we revise here elements of the standard similarity argument for plasma physics, but in a form which is specifically useful for our case.
The basic assumption put forward here is that transport physics in the considered plasma can be formulated in terms of a Fick-law expression that relates the net heat flux Q to the resulting temperature T profile: where C 0 is a constant, λ is the local normalized temperature gradient (λ = −R∂log(T)/∂r), γ is a number which indicates the type of turbulence spatial correlation regime: γ = 1 is Bohm-like (long correlation, global), γ = 2 is gyro-Bohm-like (short correlation, local). Of course one could fix γ = 2 a priori, and let the deformation from gyro-Bohm to Bohm scaling come from ρ glob * .Ξ is the dimensionless effective energy diffusivity. Note that we also assume Larmor-scale turbulence to dominate over collisional transport (neoclassical transport). ρ glob * indicates finite-Larmor-radius effects in turbulence development connected to global effects [7][8][9]. Now, to have self-similar plasmas, we need to fixΞ, that is all its parameters. This means, we fix the values of {ρ * , β, ν * , q, λ}, where ρ * = ρ s /R is the ratio of the sound Larmor radius ρ s = √ mT/(eB) to the system size R, β = 2µ 0 nT/B 2 is the electromagnetic parameter, ν * ∝ qRn/T 2 is the collisionality parameter and q is the safety factor. Moreover, the plasma boundary shape parameters (elongation, triangularity) have to be also matched. Finally, one equates the heat flux to the input power P = Q. Translated to plasma and engineering dimensional parameters, this means fixing the values of {I p , P, B, n, T} = f(R), with I p the plasma current, B the magnetic field, n the plasma density, as a function of R the plasma size.
However, following the similarity argument just introduced, one gets that, for example, P ∼ R −3/4 , thus it is not possible to scale from a bigger machine to a smaller machine, since the power would increase prohibitively. Moreover, there are some limitations in the physics description as follows: The limitations are linked to the plasma dynamics itself, since the similarity argument deals only with a specific plasma, assumed in stationary conditions. Transient phenomena are not usually treated in that context (although they are derived from the generic continuity equation with time-dependence). Additional regime-relevant parameters like, e.g. the temperature ratio T e /T i , velocity shearing rates and others may also be included as well but are not typically part of the standard argument and its applications.
Here, a different approach is suggested, which aims at obtaining a scenario, which is as close as possible to the one expected in a future reactor-relevant plasma. First, one sacrifices those parameters which, via a series of arguments, are assumed to have a more marginal role in dictating the dynamics and the plasma phenomenology inside some ranges that shall be defined as well. Second, a way how to implement practically the new approach will be presented.

A modified similarity argument: new scaling parameters
Here are defined the parameters that are chosen to scale from a reactor scenario. First, starting from formula (1), we drop the dependence of the normalized transport coefficient Ξ on the parameters ρ * , ν * , β, q, which rationale will be given later. As such, one gets: we now assume the plasma turbulence structures to be in the gyro-Bohm local regime (γ = 2), and we choose to fix the value of λ (profile peaking), and thus also ofΞ. The meaning of this latter operation is the following: we assume that between the reactor scenario and the present experiment, the profile are self-similar. That is, they have the same normalized inverse length scale λ. Since we dropped all the other dependencies ofΞ except for λ, this means that normalized transport is also self-similar, i.e. the same. As such, we can drop it from the scaling argument. Upon substituting the formula for ρ * , we arrive at the following constancy relation: with fixedQ, which can be understood as a 'normalized' heating. Again, we now equate the input power to the heat flux P = Q. This leads to a link between (P, T) as a function of (n, B). Here it is clear also while we choose the gyro-Bohmtype of scaling parameter γ: not only it is justified by assuming that turbulence spatial correlation assumes the local character by going to bigger machines (in practice, a useful rule of thumb is ρ * < 1/200 [8]), but also leads to the disappearance of the major radius from the formula ofQ, which is an advantage in terms of matching the value with a smaller machine. However, one needs an additional scaling to separate the dependences of (P, T), sinceQ = const leads to a curve P ∼ T 5/2 ; instead, we want to pin down a specific couple (P, T) ref which can be used as a reference. As such, as a second criterion/scaling formula, several choices can be made: for example, one could employ a scaling for the H-mode pedestal pressure, assuming H-mode operation. In this case, the parameters in equation (3) are intended as evaluated at the pedestal top. However, in our case, we make a more global assumption which is empirically driven. In the latter case, the parameters appearing in equation (3) are considered as averages over the plasma radius.
Finally, the second chosen parameter is: which can be understood as a measure of the average plasma transport. Fixingχ is equivalent to assume assuming the confinement time to scale dominantly with the plasma crosssection ∼ R 2 . The rationale(s) behind the choice ofχ as scale parameter are twofold: on one hand, we follow the step-ladder to a reactor which goes to the DEMO tokamak reactor [10], of size R ≈ 9 m. Plugging in typical DEMO parameters [11,12] in the ITER scaling law [13], and comparing to typical AUG parameters, leads empirically to a scaling τ E ∼ R 2 . As such, this scaling is not based on first principles, but rather is seen as a consequence if one believes that the ITER98(y,2) scaling captures the confinement trend from AUG to DEMO. Note that according to the standard identity argument, the confinement time should scale as τ E ∼ 1/B for perfectly matched plasmas. On the other hand, a constant average transport is also a coincidental observation of the L-H transition power scaling, assuming that the L-H transition happens at a constant value of the edge E × B velocity [14,15]. The logic chain is: The observation that χ ≈ const is so far of empirical nature rather than theoretical, although the derivation can follow physical intuition as done in [16].
Notice that neitherQ, norχ are expressed taking into account numerical constants or unit conversion factors. This is because what matters is that their values are kept constant between the reactor and the present machine, not so much the value per se. As such, only the actually varying variables are evidenced.
Finally, upon equatingQ andχ between the hypothetical reactor values, and the present experiment, one arrives at the following scalings for the reference applied power P ref and the reference temperature T ref as a function of the size R, the density n, and the magnetic field B, as: with C P =Q(χ/Q) (5/3) and C T = (χ/Q) (2/3) . As expected from the previous discussion relatingχ ∼ const to a similar observation done at the L-H transition, the scalings found are closer to the classic Martin scaling for the L-H transition power for P [17]. The meaning of these two scalings (5) is the following: if the plasma obtained in the present experiment, at injected power P = P ref , has, for example, an average temperature T below the reference value obtained with the scaling, i.e. T < T ref , this means that the projected scenario will have a lower temperature as well (that is, lower output fusion power). Obviously this argument is more qualitative than quantitative at this stage, taking into consideration the many reductions and approximations. In particular, we discuss now the role of collisionality ν * and electromagnetic parameter β that we have neglected.
It is well known that ν * affects plasma confinement in many ways. Above all, it has a strong effect on the density profile [18,19], and a direct effect on the electron branch of electrostatic micro instabilities, the trapped electron mode (TEM) [20,21]. It also affects the bootstrap current [22] and the poloidal flow [23]. Nevertheless, we know that all these effects saturate for ν * → 0.
Regarding β, the issue is somewhat more complicated. Electromagnetic effects influence both microturbulence [24][25][26][27], and magneto-hydro-dynamic (MHD) [28][29][30][31]. Moreover, we usually operate in a range of values which can be lower or in proximity of the reactor regime. Since at this stage we try to limit the complexity of the problem at hand, for the moment this parameter is ignored in the following. In a future work, it is planned to systematically compare the ranges of both ν * and β that are achievable simultaneously in AUG, with the perspective of using candidate discharges as reference for the reactor emulation experiments.
In conclusion, in our opinion this is a reasonable approach to drive an experimental realization of such a scenario. Moreover, one could experimentally push the non-matched parameters as much as possible using the available resources, or perform a scan around the available parameters and obtain an indication on the trends. This systematic work is planned for the future both in terms of discharge data mining and proposal for new experiments. Another aspect is related to the emulation of plasma transients as well, and the non-linearities arising in burning plasma scenarios where α power is the dominant heating channel. The proposed scalings allow to do this, by setting two main reference values (P ref , T ref ), around which we can emulate the main characteristics of the reactor plasma, as detailed in the next sections.
Finally, coming back to parameters that we wish to match: • safety factor q and plasma shaping parameters; • ratio of heat flux into electrons to the total one Q e /Q; • Ratio of particle source to confinement time; • Negligible torque input.
Matching of q is not considered as a priority (although in AUG typical values between ≈4. and 4.5 are obtained in standard H-modes). More focus is put on the second parameter, that is electron heating fraction, which is important to determine the ratio of electron to ion temperature. Here a bit of discussion is worth on why we prefer to control Q e /Q over T e /T i . First of all, since we wish to emulate dynamical phases and transients, we want to let the temperature evolve consistently with the input fluxes, and not vice versa. The idea is that the outcome in the temperature ratio should then reflect the input flux ratio, and consistently lead to a transport regime similar to the one at play in the reactor plasma. Second, the parameter Q e /Q i (so that Q = Q e + Q i ) has been shown to be a good proxy for transport characteristics [32]. However, from a control point of view, we consider Q e /Q, which value has the range [0, 1], better suited. Note that by setting Q e /Q = x, one gets Q e /Q i = x/(1 − x), which means that Q e /Q i is monotonically proportional to Q e /Q. Lastly, particle source and torque input are expected to be low in future big-sized reactors. The main reasons are the following: with respect to fueling, neutrals coming from the edge are not expected to penetrate far into the plasma, due to the high opacity of the same. Moreover, pellet deposition profiles are expected to be peripheral; finally, a pure-ECRH system is considered. As such, the core region can be considered source-free in terms of particles. Regarding torque, mainly the absence of neutral beam driven torque is the reason to consider it small. However, since in the following we will adopt neutral beam injection (NBI) to emulate most of the power sources, it is not obvious that both particle source and torque would indeed be negligible. On the contrary, they could have an impact on the various profiles, especially the density profile [33].

A fully empirical scaling set.
Before moving to the emulation of the different power sources, we add a concrete alternative method to scale down the reactor scenario. We propose to use the constancy of confinement factor H = τ E /τ scal , with τ scal = τ ITER98(y,2) , and of the ratio of power to L-H transition power: Q/Q LH . This means that, at constant plasma shape and safety factor q, we can write: and I p ∝ RB. After some manipulations we get the new set of scalings: which, aside from slightly different exponents and a weak density dependence in T ref , look very similar to (5). Interestingly, although the ITER98(y,2) scaling law has a strong current dependence but a weak field dependence, a characteristics shared between H and L modes, where for the latter a fully theory based explanation has been given recently [34], we get back a full field dependence in T ref , by choosing to work at constant safety factor q.

α power emulation
In a reactor plasma, α power will be the dominant plasma heating. That is: P α ≫ P aux (400 ≫ 0-40 for DEMO and 80 > 50 for ITER). P α depends non-linearly on the plasma DT bulk temperature T i , as well as quadratically on the fuel density.
To capture these dependencies in a meaningful way, a reference temperature is needed to be able to scale down the fusion cross section. Moreover, the distribution of P α onto ions and electron species is also a function of the bulk ion temperature. We where p α,0 is a numerical constant, n is the instantaneous plasma density (evaluated as the line average density), whereas n ref is the reference density of the selected plasma, which is the same value that is used to evaluate the reference power in equation (5). Moreover r e () is the DT fusion reactivity, and T fus is the temperature of the fusion reactor plasma. Specifically assuming the DEMO profiles given in [35], one gets T fus /(T ref ) r ≈ 2. When applied to AUG then we can write (and n = n ref ): such that for T = (T ref ) a , we get P α = (P ref ) a = p α,0 (with P aux = 0). This is how to determine the value of p α,0 . Later a concrete example with numbers for a typical AUG standard H-mode will be given.

Auxiliary, equipartition, and radiated power emulation
Auxiliary power is expected to be used in a reactor as both a mean to achieve the H-mode operation, but also, during the burning phase, to control the fusion power level or to counteract unexpected transients. The foreseen 'continuous' power from auxiliary sources is roughly 1/10 of the foreseen α power (∼40 MW as opposed to ∼400 MW of α power). On the other hand, the power needed to enter into H-mode is estimated to be about 2-2.5 times larger (up to ∼100 MW, neglecting selfconsistently generated α power). These ratios are kept for the scaled down experiments as well. Notice the important fact that since the self-generated power P α is dependent on the plasma temperature via the reactivity r e (T), it is not possible to arbitrarily choose P ref and T ref to make the plasma of arbitrarily high performance. In this sense, the fact that α power has a degraded dependence on T i at increasing temperatures, has a 'beneficial' effect on making the scaling argument sensible.
Regarding the equipartition power (i.e. the collisional exchange of energy between species), we note that we want to be able to control the heat flux in each channel (electrons and ions), such that we can emulate the actual shared power in the reactor. The ratio of equipartition time to energy confinement time is known to scale favorably for coupling, towards the reactor [35,36]. That is, correctly estimating the equipartition power is fundamental to be able to predict the ratio of the electron to ion temperature τ = T e /T i . In the following we propose a possible way to scale this down from the expected value in a reactor, based on physical arguments.
First of all, we express the classical equipartition power as (ignoring impurities): where ln(Λ) is Coulomb logarithm, Z, A are the ion species charge and mass in AMU, n e the electron density in 10 19 m −3 , temperatures are in keV, and the power is in MW m −3 . C 0 is a numerical constant that is defined below. This power is customarily subtracted from electron power source and added to the ion power source. The constant C 0 is determined such that, for a typical Hmode scenario, the ratio of collisional (τ coll ) to energy confinement time (τ E ) is scaled to be the same as in DEMO, that is: For 'DEMO = AUG', C 0 = 1. Since (τ coll /τ E ) AUG / (τ coll /τ E ) DEMO > 1, emulated equipartition is stronger than the real one. As shown in (11), the (τ coll /τ E ) quantity parametric dependencies are: (τ coll /τ E ) ∼ T 1/2 P/(n 2 R 3 ). In the case of the reactor, the cubic radial dependence is the decisive factor that makes equipartition power growing more efficient.
Finally, the radiated power in a reactor is supposed to have a dual role [35,37,38]: reducing the power crossing the separatrix, and protecting the divertor by achieving a detached regime. For our purpose, we consider only the former radiation contribution, and suppose that the radiation fraction in the core plasma is the same in both the reactor and the present experiment. Clearly, the issue of heat exhaust outside of the plasma core is of high importance for a reactor. In practice, the typical recipe for an H-mode is trying to operate the core plasma as close as possible to the L-H transition loss power, whereas an additional seeding impurity is puffed in the divertor region, e.g. argon, to provide substantial radiation losses and achievement of the detached regime [39][40][41]. In a future work, this new element will be added to our simulation tool.

A practical example using AUG parameters
In this subsection we demonstrate how to apply the proposed scalings. First, we define the hypothetical reactor parameters, which we identified with the pedix 'r': Moreover (P aux ) a ≈ 0.8 MW and (P α ) a ≈ 6.2 MW; the latter value is the one to be used for p α,0 in formula (9).
The latter has the following meaning: in the case in which n = 6 and T = 1.5 keV, as per the references, then the emulated P α = 6.2 MW as prescribed. Variations around this value are properly taken into account by the α power expression and the cross-section dependence on the temperature.

Duality of heating sources
To conclude the REP recipe, one now needs to assign the scaled down power fluxes to the actual heating system of the present experiment. We assume here that electron cyclotron resonant heating ECRH (EC) and NBI are used as heating methods (plus Ohmic power OH). Plasma impurity radiation is also considered (which can be different between the present experiment and the reactor radiation scaled down). Below, the dual expressions are presented, which are also one of the key new results of this work: (P e1/2 ) a = P EC1/2 + P NB,e1/2 + P OH1/2 − c r P rad − P eq (P T1/2 ) a = P EC1/2 + P NB1/2 + P OH1/2 − c r P rad (P sep ) a = P EC,abs + P NB,abs + P OH − P rad (P e1/2 ) r = P α,e + P aux − c r P rad,scal − P eq,scal (P T1/2 ) r = P α + P aux − c r P rad,scal (P sep ) r = P α + P aux − P rad,scal (12) where the label '1/2' means: absorbed power integrated up to mid-radius. While 'abs' means 'absorbed' in the entire plasma.
The Ohmic power at mid-radius can be evaluated as half of the total Ohmic power, or entirely neglected if only the burning phase is of interest. The factor c r represents the fraction of radiated power inside mid-radius. This can be set as a constant number (as such it is an approximation of reality). With these definitions, we make the duality concrete, by imposing the following operational equivalence: That is, we force the electron-to-total heating ratio at mid radius to be the same between present experiment and reactor, and then we force the net power into the system to be the same as well, where the 'reactor' one is scaled down using the expressions shown in the previous sections. Finally, we conclude by establishing the way how to enforce relations (13): by means of actuation. In practice, a MIMO controller is employed to steer the EC and NBI powers such as to achieve those two goal values; i.e. two actuators for two control variables.

A remark on heating power profiles.
Since ECRH+NBI will not match the α+aux powers expected in the reactor, profile-wise, this is why we prefer to match the ratio of electron to total power, and not just the electron power alone, since the electron power alone can be very different in amplitude compared to the one expected in the reactor. The drawback of this side-effect is that the more peaked ECH profile in particular will influence the plasma profile gradients and thus confinement itself.
A way out of this issue would be to perfectly match the heating profiles as well. While NBI is naturally broader than ECH, thus getting closer to what one would expect for the α power deposition profile, for ECH a proper solution needs to be addressed. One possibility could be to use multiple gyrotrons, aimed at different radial positions (e.g. 0.1, 0.2, 0.3 for example), each one modulated such as to provide the incremental integrated power that would look the same as the one expected in DEMO. This of course only for the α power fraction. Since in DEMO ECH is expected to be used anyway, that contribution to the total power would still be extremely localized in the plasma core. The dedicated solution to the heating profile problem is left for future work.

Application to an AUG discharge using Fenix
In this section, we show how to practically implement the protocol proposed in section 2, using AUG as the experimental device in which to carry out these experiments. Note that a few things need to be adjusted to be compatible with the machine capabilities, in particular regarding the available information from real-time diagnostics. Specifically, for the estimation of the ion temperature needed for both the α power and the equipartition power, we adopt the following recipe. Upon collecting the real-time values of the total plasma energy W and the electron component W e (obtained combining density and ECE real-time information), we estimate the ion component simply as W i ≈ W − W e − W fast , where the latter contribution W fast is computed by the real-time code RABBIT [42]. Then, the equipartition power is replaced by this approximated expression: where C 1 is a numerical constant. This constant is chosen as to fit the evolution of the actual equipartition power on average.
The real-time information on the heating profiles is obtained via both RABBIT and TORBEAM-rt [43,44]. The radiated power is obtained from selected bolometry channels [45] which should return a real-time estimate of the core component of the total emitted radiation.

Experimental scenario
To test the proposed protocol, standard H-mode discharge #40446 is adopted as the prototypical 'burning plasma' reference. This discharge parameters are: plasma current I p = 1. MA, magnetic field B T = 2.5 T, average density ⟨n⟩ ≈ 6 − 8 × 10 19 m −3 . The actual discharge evolution is displayed in figure 1. It can be noticed that the chosen discharge operates at a rather high density. To carry out the REP simulation later, we artificially reduce the density by a factor of ≈1.5. The reason for this choice is the following: standard H-mode discharges at AUG are performed at rather high-density, which does not allow enough flexibility in variation of the Q e /Q parameter. By artificially reducing the density, we allow for more range of variation in which it can be controlled. In any case, when doing planned experiments for this specific emulation protocol, it is better to work starting from a medium-to-low density plasmas.
To study the REP protocol proposed in section 2, the experimentally used heating scheme is fully replaced by the recipe presented above, whereas, in the ASTRA code [46,47], a simple controller of the heat fluxes is implemented, to achieve the needed power levels given by (13). In particular, NBI is used to control the total needed power Q, whereas ECH is used to control the electron power fraction Q e /Q at mid-radius. That is, the control protocol is the following: with K T P , K T I , K E P , K E I are the four proportional-integral control constants (positive numbers adapted to get the best response). Finally, the NBI power is obtained simply as: P NBI = P T − P EC . Obviously both EC and NB powers are capped at 0 from below, and the maximum available power is limited by the actual actuators capabilities. For the simulation shown later, we fix the PI constants to these values: K T P = 1, K T I = 50, K E P = 25, K E I = 250. Regarding the actuation itself, we employ a realistic description of the NBI and ECRH actuators, assuming a modulation duty cycle of 20 ms, whereas the power is delivered in a shorter window, such that the average injected power matches the requested power. The rise time of the heating beams are respectively 3 ms for NB and 1 ms for EC. As a related note, the main reason why it is better to use a control scheme instead of assigning directly the requested power to the individual sources, is exactly because of the many subtle sequences from 'beam on' → 'power absorbed in plasma'. The flexibility of the PI controllers is enough to allow this better approach.

Set-up in Fenix
The tokamak flight simulator Fenix [48][49][50][51] is the optimal tool to carry on these investigations. Fenix simulates the entire discharge from coils pre-magnetization to plasma termination. It consists of a coupling between a virtual copy of the ASDEX Upgrade Discharge Control System in Matlab-Simulink TM and the ASTRA transport code. This coupling allows to perform full-discharge simulations including the response of the control system, and of the plasma physics from the side of the ASTRA transport solver. It includes reduced models for all the plasma physics aspects, at least at the most basic level, but considering known non-linearities in e.g. fueling, heating, and core transport.
As said previously, for this test, we employ discharge #40446, but instead of the power trajectories prescribed from the pulse schedule used in the actual experiment, we replace them with the control scheme of (15), where the reference signals are given from the duality formulas (12) and (13). For the example shown here, we adopt the following choices for some of those quantities: the equipartition power scaled is set as P eq,scal = 10P eq , with P eq defined by expression (14). The value of 10 is chosen to replicate the stronger species temperature coupling expected in a bigger device (the scaling factor can vary between ≈3 and 10 typically). Regarding radiation, we set P rad,scal = P α + P aux − 4, so that (P sep ) r = 4. MW. This means that the heating level in the simulation is such as to add up to 4 MW (subtracting the actual radiation given by the background plasma). For comparison, in the reactor, P α ≈ 400 MW, P aux ≈ 0, and the P sep has to be close to the L-H transition, that is P sep ≈ 150 MW. As such, P rad,scal (reactor) ≈ 250 MW. As discussed in section 2.4, this radiated power would be produced by seeded impurities (e.g. Xe, Ar). For the present experiment (virtual experiment in Fenix), we plan to replace the radiated formula with P rad,scal = P rad (the latter being produced by the actual background plasma), and instead injecting impurities such as to get (P sep ) r = 4. MW. However this is not done here since a model for impurity seeding is not yet available in Fenix.
Regarding the formula for the α power, equation (9), we rewrite it as: with p α,0 = 6 MW, and W being the diagnosed plasma energy, The control scheme is active all along the discharge. To push the plasma into H-mode, P aux = 3 MW is set for t = 2 to t = 3 s. Otherwise it is set to 0 (assuming 'ignited' scenario). Moreover, we neglect the density dependence in the fusion power formula (8), since in the simulated discharge density is not controlled and has a few large excursion as shown in figure 1.
To carry out these simulations, the employed core transport model to predict the temperatures is a simple gyro-Bohm model, calibrated to give a confinement time roughly of the same size of the ITER(98,y2) scaling time.

Results and discussion
The results of the simulated 'burning plasma' scenario are given in figure 2. In the figure, we show time traces of several quantities of interest. In particular, in the top-left panel, the actuator powers P EC , P NB , which are controlled based on the request powers (15), are compared to the calculated separatrix power loss P sep , and its target P sep,t . On the top-right panel, the electron heating fraction is compared between the target request 'target', and what is actually obtained in the modeling 'sim'. In the bottom-left panel, the modeled electron and ion temperatures are shown, the value taken at the magnetic axis. Finally, in the bottom-right panel, we show the calculated α power, which is used as part of the total power request as per (12), and its sub-division in electron and ion. Moreover, the 'auxiliary' request power and the 'radiated scaled' power are also shown. These also enter into the total reference power.
The REP is activated from t = 2 s forward. It can be seen that both the target P sep and the target Q e /Q are followed nicely by the simulated values. Since the partition of powers is such as to get that electron heat flux ratio, the temperature applied EC power, NB power, measured separatrix power, and target separatrix power. Top-right: target electron heat flux ratio at mid-radius, and simulated one (red). Bottom-left: electron and ion temperatures on axis. Bottom-right: calculated α power, electron part, ion part, auxiliary power, and requested radiated power. ratio comes out self-consistently as can be seen in the bottomleft panel. Finally, the calculated 'reactor scaled down' values of α, auxiliary, and radiated power are shown. Notice that we are not injecting any impurity (as opposite to what would be actually done in DEMO) to manage the requested radiated power. Instead, this is directly subtracted from the requested EC power. In this way, we avoid the need to use an impurity transport model. Of course, the limitation of this approach is that all the time scales and profile effects of the impurity are ignored.
For summary, we collect the equivalence between the 'reactor' and AUG parameters used in this simulation in table 1.
Note that the ECRH power depicted in figure 2, which is observed to be about 2 MW at maximum, will increase if the radiation power is obtained from impurity injection into the actual plasma, instead of being fixed from the scaled radiation. In the case of seeded impurity radiation, ECRH power request can go up to 4-5 MW for this case.

Conclusions
In this work we have presented a practical protocol to emulate a fusion reactor scenario, having in mind as basis the standard H-mode type of burning plasma scenario. This is achieved by assuming two specific scalings, as well as a recipe to emulate all the various types of applied power, and power absorbed by the different species, in a reasonable way.
Using the flight simulator Fenix to carry out the virtual experiment, it is shown that this protocol works, and both the total heating and the electron fraction are reproduced, leading to a realistic prediction of the temperature ratio as a byproduct.
Future work is devoted to two aspects. From the physics point of view, we go into more details into how to best match most of the parameters that we can in AUG, such as to obtain the closest possible replica of a reactor plasma. From the operational point of view, replicating the control scheme in the control part of Fenix, in particular modeling the actuator technological capabilities in a more detailed and accurate way.