Development and application of a predictive model for advanced tokamak scenario design

Advanced tokamak (AT) scenarios applying additional heating during the current ramp (early-heating) usually require many iterations if developed fully empirically. To reduce the required experimental time, a model has been developed in the ASTRA framework, capable of doing predictive simulations of the relevant parameters. As scenario development requires fast iterations and inter-discharge runs, a sufficiently short run-time is required. While using a simplified transport model to achieve this, comparisons to experimental data from ASDEX-Upgrade (AUG) still show good agreement. Using this model, a new high performance early-heating AT scenario has been developed and successfully run on AUG with the results matching the predictions.


Introduction
The inherently pulsed operation of conventional tokamak scenarios, where the plasma current is driven mostly inductively by the central solenoid, is not desirable in future fusion power plants due to cyclic loads on the machine and for economic reasons due to down-time in between discharges.Advanced tokamak (AT) scenarios [1][2][3] offer a solution to this situation by maximizing non-inductive current-drive.This 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. is done by using external heating systems to drive current offaxis, thereby raising the safety factor (q) profile such that an increased bootstrap current (j bs ∝ q∇p) fraction is achieved.The resulting setup yields an increased maximum pulse length up to potentially fully non-inductive operation [4].A positive side effect is a potentially improved stability and confinement since low helicity magnetohydrodynamic (MHD) modes, such as the sawtooth instability, are suppressed, when their resonant flux surfaces are no longer present due to an increased q min .However, the broader current profiles due to the elevated q min potentially reduce the ideal stability and thereby the maximum achievable β N [5].
Most of the recent studies at ASDEX-Upgrade (AUG) [6,7] have used this approach.However, in such a setup the q-profile can drop considerably below the intended value before the actuators are applied.The subsequent increase of the q-profile develops with the current diffusion time τ r = µ 0 σL 2 , where σ is the volume averaged neoclassical conductivity and the characteristic length L is usually set to half of the minor radius.For a typical AUG AT scenario it is ~1.5 s, making this approach acceptable.However, in future larger machines such as ITER a much larger current diffusion time of a few hundred seconds (estimated based on values from [8]) is expected, requiring an alternative solution.
The alternative, 'early-heating' has the additional current drive already applied in the ramp-up phase of the plasmacurrent, avoiding the intermittent drop of q below the target value.As the loop voltage from the transformer coil is applied at the outer edge of the plasma, the inductive current arises there and then diffuses inward on the timescale of τ r towards a peaked current profile where j ∝ σ ∝ T 3/2 e .The early heating has the additional benefit of increasing the conductivity and thereby the current diffusion time of the early plasma, thereby slowing down the diffusion of the current towards the plasma center.The resulting plasma current retains a broader current profile, corresponding to a higher q-profile.The downside of this setup is the high sensitivity towards the timing of heating and fueling systems as well as a reduced stability towards tearing modes, which can give rise to stability issues, causing fully empirical scenario design to usually require many iterations.
Due to the difficulty in designing such a scenario, while early heating scenarios have been studied on Carbon-AUG in the past [15,16], recent studies have put less focus on them.This paper addresses this issue by providing a workflow to generate a stable path to an early-heating scenario in a full tungsten machine.
Section 2 gives a brief overview of the experimental setup for AT scenarios at AUG.In order to move away from fully experimental scenario design, a model capable of predicting the behaviour of relevant quantities has been developed.This model is presented in section 3.In order to validate this model, it was used to design a new early heating scenario for AUG.This process is shown in section 4. Performance of the new scenario is discussed in section 5.A summary of the results and an outlook is provided in section 6

Available current drive systems
AUG features two neutral beam injection (NBI) systems with four beams each, capable of producing 2.5 MW per beam for a total of 20 MW [17].The paths of the beams through the plasma are shown in figure 1.Of those, the two beams with a considerably longer path in the plasma (six and seven), which are also the furthest away from the magnetic axis (and therefore called off-axis beams) generate the highest amount of off-axis current drive.The other beams are considerably closer to the magnetic axis (called on-axis beams) and therefore less efficient at driving off-axis current, although 2 and 5 still are considerably better than the rest.The driven current  NBI driven current as calculated by RABBIT [18] in ASTRA [19,20] for each source individually, note that in the real discharge (see section 4.1) not all beams were running.per beam as calculated by a simulation is shown in figure 2 for an example case.Note that in the real discharge not all beams were active.For technical reasons, one of the off-axis beams (number seven) was not available for parts of this work.
Additionally AUG has eight gyrotrons for electron cyclotron resonance heating (ECRH) with a power of up to 1 MW each [21].The deposition location of all gyrotrons can be moved, even during a discharge by changing the angles of their injection.The same process also allows to change how much current is driven, down to effectively zero.As the current drive efficiency depends on T e /n e and T e is considerably more peaked, it is reduced when moving towards the plasma edge.

Typical AT operation at AUG
Figure 3 shows a comparison between a current profile, as it can be seen in a low-power H-mode on the top left and an advanced scenario on the top right.In order to move away from a centrally peaked current profile, which is dominated by inductively driven current, the q-profile is increased through the use of NBI and EC driven current to increase the bootstrap Figure 3. Current profile during a q 95 ∼ 5.2; Ip ∼ 800 kA advanced scenario (right) with hatched regions denoting a negative contribution compared to the current profile in an arbitrarily chosen H-mode (left) and the corresponding q-profiles (bottom).
current.The bottom figure shows the q-profiles corresponding to the two current profiles.In the advanced scenario case the majority of the externally driven current is provided by the NBI systems, while the lower ECCD fraction is important to shape the q-profile to allow for a high bootstrap fraction.It can also be seen, that in this process the maximum of the current is moved away from the axis.
In order to maximize the available external current drive, AT scenarios are run at low collisionality.This requires operation close to a boronization [22] to minimize impurity influx.Additionally, a plasma shape with high large wall clearance is chosen to minimize tungsten sputtering from the wall [7].It has been shown in [23] that the divertor neutral density (n 0,div ), which is correlated to the separatrix density has a significant impact on global confinement.Here the global confinement scales inversely with n 0,div , which is therefore minimized.For that reason density control is done via feedback on this quantity.
All scenarios in this work were run with a feed-back controlled flat-top plasma current of ~0.8 MA, yielding q 95 ∼ 5.2 and a magnetic field of ~2.5 T. These values were chosen to be comparable to estimations for an ITER steady-state scenario [24].
In devices with tungsten plasma facing components, sufficient fuelling has to be provided before the L-H transition in order to minimize the ELM-free phase after the transition, which would otherwise lead to problematic tungsten accumulation in the plasma core [25].This puts a lower limit on the plasma fuelling.An upper limit to the density is created by the possibility of the appearance of a MARFE [26], if the density gets too high before sufficient heating is applied.This yields a rather limited density range, where the plasma reaches a stable H-mode.
In order to prevent tungsten accumulation in the core, one gyrotron is used for central heating [27] with no current drive.For sufficient current drive to reach relevant q-profiles, one of the off-axis NBI systems is required at all times.The other NBI sources can add additional current drive and increase β, however due to the low density, they can not always operate for the entire discharge duration, due to overheating of the target tiles on the inner wall from shine-through.The remaining gyrotrons are used to do adjustments to the q-profile or help with NTM suppression [28].

Relevant diagnostics
Electron temperature (T e ) and density (n e ) profiles are reconstructed by the integrated data analysis (IDA) [29] tool, combining data from the electron cyclotron emission (ECE) [30], Thomson scattering [31], DCN Interferometry and the Lithium-ion beam emission spectroscopy [32] in a Bayesian approach.
The effective charge is calculated based on line integrated bremsstrahlung measurements from CXRS together with temperature and density data from IDA.Since only a few lines of sight in the core can be used, due to issues caused by reflections, usually only an averaged value is provided [34].
The integrated data analysis equilibrium (IDE) [35] tool generates a plasma equilibrium by solving the Grad-Schafranov equation based on kinetic electron and ion profiles, fast ion density, and rotation profiles constrained by the pressure profile, data from magnetic measurements and current diffusion.Profiles are converted to a base of normalized toroidal flux (ρ tor ) using this equilibrium.
It has been shown in [36], that current diffusion fits well to internal measurements of B pol .If these internal measurements are not included, the error bar gets considerably larger, but the q-profile does not change systematically.For the experiments presented during this work, internal current measurements were not available.Based on the previous observations, we expect current diffusion to describe the q-profiles sufficiently well.This statement is supported by the fact, that MHD markers, found during this work are compatible with the q-profiles and no large MHD activity, that would lead to a redistribution of current is present.

Model
In order to guide the experimental scenario design, a model has been developed in the ASTRA [19,20] framework, which can simulate the behaviour of the relevant parameters.As a 1.5D transport code ASTRA calculates a 2D equilibrium and a 1D diffusion equation.A schematic overview of the model is shown in figure 4: Initial input is the actuator setup.Since AUG has only ECRH and NBI heating (see section 2.1) as heating systems capable of supplying relevant amounts of current drive, so far only those two systems are considered.The modular design of ASTRA does allow adding additional systems in the future.The respective heating and current drive is calculated by the TORBEAM [37] code for ECRH and the RABBIT [18] code for NBI.Both of these codes are coupled in ASTRA such that they use the output from the previous time-step as input.
The next input is the plasma density.Since it is feedback controlled (see section 2.2) in all discharges considered during this work and does not change significantly between reasonably similar discharges, experimental data from reference discharges are used here.In order to have a sufficiently high density at the L-H transition, the time axis of the density is changed such that the density at the time, when the first NBI system is turned on is the same.
A heat transport model is required for the calculation and has to be provided.The exact model used is described in section 3.1.ASTRA then uses these inputs to calculate temperature profiles.Based on these and the current drive input the current diffusion equation is solved, yielding a current density profile, which is used to calculate the q-profile.

Transport model
Since the intended use for the model is scenario development, where iterative testing of various different setups is foreseen, and running the model in between discharges may be of interest, a short run-time is required.Note that most of rampup (starting from t = 0.15 s corresponding to ~0.3 MA) plus multiple current diffusion times of the stationary phase (usually up to t = 4.5 s) are modelled.
While there are first principle models (such as TGLF [38] or QuaLiKiz [39]) available, using them would have a run time on the order of hours or longer.For the purpose of this work, a faster model is required.Therefore, a comparatively simple, fully analytical Bohm/gyro-Bohm model is used, achieving a run time of a few minutes on contemporary hardware (Intel Xeon Gold 6130, using 8 cores).This model is adapted from [40][41][42].Multiple free parameters c n are present in the equations.The Bohm/gyro-Bohm electron heat conductivity (χ e ) is: which includes a simplified TEM threshold: The Bohm/gyro-Bohm ion heat conductivity is: which includes a simplified ITG threshold [43]: with the magnetic shear S and the effective ion charge Z eff .It has been found empirically, that including a fast ion term f FI and an electromagnetic term based on zero order results from gyrokinetic simulations [44] improves the agreement in the ion temperature: with the fast ion pressure P FI , the toroidal magnetic field B tor and the normalized magnetic shear S norm = 1 − S max(S) .Transport in the edge/pedestal region of ρ tor ∼ 0.95-1.0, is included through a scaling law, which calculates the electron heat conductivity at the pedestal as: where P kin,ped = k B (n e,ped T e,ped + n i,ped T i,ped ) with the Boltzmann constant k B is the kinetic pressure at the pedestal and the critical pedestal pressure is given by: With the poloidal magnetic field at the pedestal B pol,ped and the critically stable poloidal beta at the pedestal β pp , given by the scaling [45]: With the elongation κ, the triangularity δ, poloidal beta β pol , the Greenwald density n gw and a machine dependant scaling parameter for the pedestal width w p , which is set to 0.11 for AUG.
The ion heat conductivity at the edge (ρ tor ≳ 0.95) is assumed to be identical to the electron heat conductivity.The neoclassical contribution to the conductivity is calculated in ASTRA, but it is comparatively small for the cases studied here.
An example for a resulting electron heat conductivity profile based on these contributions is shown in figure 5(a).It can be seen, that the gyro-Bohm term is dominating, whereas the neoclassical contribution only has a relevant contribution inwards of ρ tor ∼ 0.1.Towards the edge, an exponential-like behaviour (shown by the blue dotted line), can be seen in the gyro-Bohm term.This would lead to a sharp drop-off in conductivity, when switching to the edge model, creating a pronounced kink in the temperature profiles as can be seen in figure 5(b).In order to avoid that, and improve the agreement with experimental behaviour, a less steep transition would be preferred.Here a linear connection between an empirically determined cutoff-point and the edge is used.Empirically the best agreement was found for ρ tor ∼ 0.8.Note that this value appears to be moving to slightly lower values during a discharge.The exact dependencies are unclear and might be explored in future work.The resulting heat-conductivity term reasonably reproduces the behaviour calculated by ASTRA through power balance with experimental T e and heating profiles as input (shown with the black dashed line).
The effective charge Z eff is important, as it affects the plasma resistivity and the NBI current drive in addition to the ion heat conductivity (see equation ( 4)).It is assumed to be (a) Contributions to the heat conductivity; the blue dashed line shows the behaviour of the gyro-Bohm term towards the edge, a linear approximation towards the edge is used to improve agreement with experimental data.The black dashed line shows the heat conductivity calculated by ASTRA, when using experimental Te-profiles and heating power as input.(b) Effect of the linear approximation in the heat conductivity; the offset in the core between the case with and without the approximation could be mitigated by changing the free parameters; the insert shows a focus at the kink at the transition between edge and core model.radially and temporally constant at a value reasonably close to the experimental behaviour (usually 1.3 ≲ Z eff ≲ 2.0).The dependence on machine conditions makes it necessary to set this value manually for each discharge.Figure 6 shows the impact on Z eff on the evolution of q at three radial locations, by comparing the effect of values on the opposite ends of the range.It can be seen, that the effect is rather small and the difference decreases over time.
Since the goal of this model is to simulate an entire discharge, including the ramp-up, describing the L-H transition is required.Here it is assumed that the transition happens, once the heating power at the seperatrix P sep surpasses a threshold (P L−H ).The threshold is adapted from [46], with simplified exponents and the prefactor turned into a free parameter to better match the experimental behaviour: (10) with the electron density at ρ tor = 0.9, corresponding to the pedestal top location in H-mode and the surface area of the plasma S lat .To reproduce the experimental behaviour, the free parameter c has been adjusted such that the L-H transition does not occur, when only one beam is active (region from ~0.65 s-0.85 s in figure 7) and does occur, after the second beam is turned on (at ~0.85 s).This has been tested for all discharges  studied during this work and it was found that c = 0.0075 produces a satisfactory behaviour.Note, that for the experiments of interest, the heating power at the seperatrix increases significantly at the time of the L-H transition (after the second beam is turned on), meaning that this parameter could be larger (moving the blue curve up) and still produce a similar result.For some of the free parameters in the gyro-Bohm model a different value is used in L-mode.Since the edge scaling is not valid in L-mode, the heat conductivity at the edge is set to an arbitrary large value (χ e/i,edge = 20) and the linear approximation is not used (see figure 8).
Experimental data is used for the initial conditions, since their impact was found to be low, due to the short energy confinement time and even more importantly current diffusion time in the early phase.This finding is consistent with findings from [47].It was found that the setup is mostly independent from the initial temperature (ASTRA requires the initial profiles to be nonzero), as any difference disappears after only a few iterations.The initial safety factor is slightly more important, as its effect can be seen for a longer time, but also disappears before the heating systems are turned on.
In order to choose the free parameters, starting from an arbitrary initial guess, a scan over them was done in order to minimize the metric: where f exp is the experimental value, f m is the value calculated by the model, σ 2 is the respective variance for the experiment and the model, i is the point of interest and N i is the total amount of considered time points.A value around one is considered a good fit, a value much larger than one indicates that the model reproduces the observation poorly.The experimental variance is determined from the error ranges of the values.As ASTRA simulations are deterministic, an error estimate can not be given and is therefore set to zero.
Looking at the T e and T i time evolution at multiple radial points for a set of 800 kA reference discharges, the optimal value for each of the parameters was determined.The radial locations of ρ tor = 1/3 and ρ tor = 2/3 were chosen, since the range closer to the edge is dominated by the scaling law and the boundary conditions and the area further in the core is of lower importance due to the low volume and higher error bars in the measurements.
As the model is intended to predict the time evolution of the plasma, this setup focuses to match the experiment over time.By doing this, the impact of fluctuations and noise, as well as the impact from using real density data as input is minimized.Due to the uncertainties and the experimental noise, there is a range, in which the quality of the fit does not change (χ 2 changes by less than 10% when the free parameters are varied in this range).Averaging over the set of reference discharges, the minimum value of χ 2 achieved for the electron temperature is 1.2, which is a quite good fit.This is important, since T e has a large impact on the conductivity and therefore the qprofile, the prediction of which is the main goal of the model.For the ion temperature, the best fit achieves χ 2 = 3.6, which is considerably worse.As the impact of T i on the other parameters is low, this was deemed to be acceptable.In section 4.2, a parameter scan over T i is done in order to validate this statement.
The best fit for the reference discharges is achieved at slightly different values for each of the parameters.Table 1 shows the mean of the values, generating the best fit for each of the reference discharges.Since the parameters c4 and c5 have very close to the same effect on the profiles, c5 was set to a fixed value and the fit was done for c4.Also shown is the average range in which the parameters can be altered without significantly changing the quality of the fit (difference in χ 2 of less than 10%).
It was found, that the optimal values for the free parameters are close to each other for all of the discharge with an overlap of the region, where the quality of the fit is similar, which gives confidence towards the predictive use of the model.The last row shows the parameters, which are finally used.These differ from the mean in some cases, where it was found that the used value improves the agreement for a majority of the discharges.
It was found, that the discharges performed during this work are well described by these parameters.

Model limitations
In its current form, the model does not include experimentally found density limitations around the L-H transition (mentioned in section 2.2).This means, that a given density setup, producing a valid result in the simulation may not be achievable in reality.The density may therefore need to be tuned experimentally to achieve a successful discharge.
The linear connection between Bohm/gyro-Bohm transport and edge/pedestal transport does improve the agreement with the experimental behaviour, but does not have a solid physics basis.Considering that the cut-off value, which would achieve the best agreement can move with time, some effects may be missed here.
MHD effects are not included in the model: a simulated behaviour may not be sufficiently NTM-stable to be experimentally achievable.This situation can be improved by avoiding low-shear regions around rational flux surfaces as it has been done in [48].Further, if MHD activity has a significant effect on a discharge, such as a drop in core temperature when a large NTM is present, these effects cannot be captured by the model.However, experimental data with NTMs present were anyway not used in this study.
The transport model only includes TEM and ITG turbulence.If other phenomena have a relevant effect on the plasma behaviour, the resulting scenario cannot be reproduced.
Experimental fluctuations, especially those caused by ELMs, are not included in the model, which causes the agreement of profiles between simulation and experiment to be worse in some time points.If experimental density data is used as input, these fluctuations do propagate and have some impact on the other quantities.
Since only AUG discharges were considered and all of them were using the same plasma current, magnetic field and plasma geometry, a possible dependence of the free parameters in the model on those quantities cannot be excluded.A more extensive physics based model, that looked at a larger set of discharges (although still only on AUG) is available and might be used in future work [49].

Application to AUG AT scenarios
In order to validate this model, it has been used to design a new scenario for AUG.Starting from a reference late-heating scenario, the goal was to generate a stable early heating scenario, which reaches a q-profile as high as reasonably possible early in the discharge without an intermittent drop.

Designing a new scenario
Starting point was a 800 kA late heating scenario (AUG discharge 37722), which was part of the set of discharges, used to determine the free parameters of the model.As a first step, the time when heating starts was set.This is limited by the time, when the X-point of the magnetic configuration is reliably established.Starting earlier is not possible, since applying external heating in the limiter phase would lead to unacceptable tungsten influx.
To account for the earlier heating, the density was adjusted, such that it is sufficiently high at the time of the L-H transition.In the simulation, this has been done, by taking the density evolution from the reference scenario and adjusting its time evolution such that the density on heating start is the same.In the experiment the feed-back controlled fuelling needs to be adjusted in order to achieve the intended behaviour.The setup works as follows: in the L-mode phase, the fuelling is based on feed back controlling the line averaged density of the line of sight through the plasma core.Shortly before the L-H transition is expected, coinciding with the time, when heating starts, the feedback control switches to the edge density.During this phase, the fuelling is increased significantly by programming a steeply increasing desired density trajectory for the duration of the L-H transition.Afterwards, the feedback control is switched to the neutral gas density in the divertor for the flattop phase, i.e. the recycling flux is controlled as described in [23].Figure 9 shows the density evolution at a representative location for the reference case and the newly designed scenario.For reference the evolution of the plasma current and the times, when heating starts are shown.
An initial current drive setup was chosen arbitrarily, based on the reference scenario.The simulation was run and the resulting behaviour of the safety factor was evaluated.
The requirements were to keep q min at the highest value still able to be sustained stationarily throughout a discharge.Further, local minima in the q-profile and its time evolution were to be avoided and jumps in the first derivative due to discrete turning on of heating systems were minimized.
In an iterative approach, changing the timings and locations of the heating systems, when the safety factor drops below an arbitrary preset value or the previously mentioned criterion is no longer fulfilled, the setup was adjusted, until a satisfactory behaviour of q was achieved.
The final setup of the heating systems is shown in figure 10(a) in comparison to the reference scenario for a discharge employing the new recipe (see also next section).
The NBI-power is increased at the maximum rate, empirically found to not cause MHD events.Additional ECRH power is added shortly after the second neutral beam is turned on, in order to control the behaviour of the q-profile in the plasma center.The gyrotron for central heating is turned on during the first step.The figure 10(b) shows a comparison of the radial distribution of the driven current for the modelled scenario after all heating systems are turned on and the reference scenario at a time point, where it is stationary.As can be seen, NBI drives more current, however a lot of it is centered close to the magnetic axis and it produces rather broad peaks.This contribution is similar to the reference case, with slightly more on axis current drive.ECRH, while driving less current, has considerably more localized peaks, allowing for a more precise shaping of the q-profile.The ECRH peaks have a very similar power, the significant reduction in current drive while moving outwards stems from the reduced efficiency off-axis.The outermost gyrotron is used for NTM suppression at the q = 2 surface and has very low impact on q.The innermost gyrotron is used for core heating (see section 2.2) in a configuration not to drive current, but generates a small negative contribution.A smoothed version of the EC driven current, preserving the surface integral is also shown.The main difference to the reference scenario here is that the deposition locations have been moved outwards and a different set of NBI sources is used, leading to a higher current drive.
The smoothing shown here is the one used in the IDE code (see section 2.3) to calculate the experimental profiles.The sharp peaks in current distribution would lead to distinct features in the q-profiles, which are not seen in the experimental data however, with the available resolution it is unclear if they could be seen, were they present.

Experimental results
The designed scenario was run at AUG and a stable discharge (39221) was achieved.In order to get to that point, some tuning of the fuelling in the early phase was required in order to fulfill experimental density limitations mentioned in section 3.2).To be precise, three discharges were required to adjust the density Figure 9. Density evolution at a representative location for a late heating discharge compared to an early heating discharge in comparison to the plasma current.The density at the time, when the heating starts is at a very similar value for both discharges.setup.At the relevant point in time (at t ∼ 0.8 s) the discharge is still rather generic, such that this tuning will likely not have to be repeated for other early-heating setups.
Figure 11 shows the evolution of the safety factor in the reference case, the predicted behaviour by the model and the experimental result of running the newly designed scenario on AUG.It can be seen that the experimental behaviour matches the simulation well, while achieving a considerably higher value compared to the reference scenario.
In the following a comparison of the simulated data with the experimental results will be shown.Figure 12 shows the behaviour of the electron temperature: Excluding fluctuations (see section 3.2), the electron temperature in the simulation matches the experimental data well, both in profile, as well as in time evolution.The fluctuations in the simulation are a result from the density input.
The safety factor shown in figure 13 matches the experimental data as well, both in L-mode and H-mode.Note that some discrepancies are to be expected due to experimental fluctuations (see q time evolution shown in figure 11).The modelled behaviour of the ion temperature (see figure 14) shows a larger discrepancy to the experiment.
The trends and profile shape are reproduced.However, the time evolution of the plasma after the L-H transition is different, indicating possibly some missing physics in the model.Changing some of the free parameters, a better agreement could be achieved, although those changes are not consistent between discharges.
A discrepancy is also visible in T e (see figure 12) in the region of 0.6 ≲ ρ tor ≲ 0.8.As there are no indications of an experimental explanation such as an internal transport barrier, this behaviour is attributed to missing physics in the model.
The ion temperature affects the q-profile mostly through the ion-electron heat exchange term affecting the electron temperature.As T e still shows good agreement even when a difference can be seen in T i , the conclusion can be drawn, that the model is able to account for this discrepancy.Therefore this discrepancy in T i is not a concern.
In order to test the direct impact of T i on the behaviour of q, simulations were carried out, where the entire time evolution of T i was fixed to the experimental data multiplied with a set of factors.In these cases T e is also fixed in order to eliminate its effects on q. Figure 15(a) shows a comparison of the time evolution of q for the experimental data, a simulation where T e and T i are fixed to the experimental values, and a simulation where T i is increased by a factor of 1.5.It can be seen, the the overall impact is rather low, which is expected as T i does not affect the conductivity but only has minor contributions to the bootstrap term and the neutral beam driven current.Considering that the deviation of the simulated T i is lower than this factor of 1.5, the discrepancy was deemed acceptable.The behaviour β shown in figure 15(b) follows a similar trend to core T i , where the initial time evolution is not fully reproduced again indicating, that some time-dependent phenomenon is likely being missed by the model.In the later parts of the discharge the simulation agrees with experimental data.

Scenario performance
The goal of AT scenarios is to use the bootstrap current to reduce the reliance of the on the inductive current driven by the central transformer.This section will deal with the impact of the increased safety factor profile on the bootstrap current and thereby the overall current-distribution.
Since the bootstrap current is also proportional to the pressure gradient, which can be adjusted through β, reactions of the scenario to increased heating power were explored.

Validation case
The current distribution for the model validation scenario (see section 4.2) is shown in figure 16.This scenario operates at plasma current 800 kA, corresponding to a q of 5.2 at a magnetic field of 2.5 T. The plasma current consists of inductive current, externally driven current (NB and EC) and bootstrap current.can be seen, that bootstrap current makes up a significant part of the overall plasma current.The average value in the stationary phase is at ~41% with a standard deviation of ~4%.While ECCD is important for q-profile shaping, its contribution to plasma current is rather small.The NBI driven current is considerably larger.Due to the switch from an on-axis to an off-axis source at ~7 s, an increase in NBI driven current from ~30% to ~40% can be seen.This increases the non-inductive fraction, achieving a value of ~90%.

Probing stability limits
In order to test the β-stability of the scenario, a repeat was done with an additional NBI power up at ~3.5 s, after beta has reached its flat-top value.
The corresponding behaviour of β N can be seen in figure 17.During the ramp, an increase in confinement (H 98 , [50]) and a considerable increase in bootstrap current can be seen.The plasma disrupts ~4.5 s after the appearance of an ideal mode at β N ∼ 3.2, indicating a limit.The identification as an ideal mode [51] is based on the constant phase of the ECE  measurements (see figure 18(a)) and the very high growth rate of the mode of γ MHD ∼ 4.2 • 10 6 s −1 , calculated based on the magnetic signal (see figure 18(b)) and [51], yielding a corresponding τ MHD ∼ 1.2 • 10 −7 s comparable to the Alfven time.
Notably, the scenario becomes transiently fully noninductive before disrupting.Since the amount of off-axis NBI power does not change, no significant change in NBI driven current is observed.
An additional repeat was done with a feed-back controlled increase of β N to 2.7 after a stationary is established at ∼3.5 s.The resulting discharge is shown in figure 19.At ~5.2 s one NBI source shuts down due to excessive heat on the inner wall target.The controller switches on an alternative source however, the intermittent loss of power causes a considerable drop in the core ion temperature (see figure 20(a)), which does not fully recover.
As can be seen in figure 20(b) comparing logarithmic temperature profiles, while the T e -profile only shows an offset between a time-point before and after the event, the shape of the T i -profile changes.The temperature in the plasma core is lower after the event, while the edge temperature stays mostly the same, indicating a change in transport behaviour.Accordingly a larger change can also be seen in the T e /T i profile in the core region.
An additional increase in β N is attempted at ~6.5 s but fails with the appearance of a quickly locking NTM, effectively terminating the discharge.This occurs at a level, where the previous discharge was still stable, indicating a reduced stability threshold due to the change in transport, leading to a change in the profiles.The remaining part of the discharge has no relevance to this work.
Once again, the impact of the increase in β can clearly be seen.In the phase before the NBI switch at ~5.2 s, a bootstrap current fraction of ~48% is achieved, the average over the entire stationary phase (from ~4 s to ~6.5 s) is at ~46%, excluding the drop.The small change in bootstrap current during and after the event shows a behaviour very similar to the confinement and T i , indicating that these changes are likely driven by them.Overall a non-inductive fraction of 87%-90% is achieved.Due to technical difficulties, the second off-axis NBI system was not available for this discharge.Assuming a similar increase in NBI driven current (of ~10%), as can be seen in the 'baseline' case (figure 16), this scenario may have become fully non-inductive with both off-axis beams active.Testing this possibility may be of interest for future campaigns.

Summary and outlook
A model was developed in the ASTRA framework, which can predict the behaviour of temperature and safety factor time evolution depending on the heating setup for an AT discharge on AUG.Currently the plasma density is a required input, but including a density model is planned for the near future.This makes it a useful tool for assisting in development of early heating scenarios, by giving information about q-profile shaping in the ramp-up phase and allowing to choose heating systems such that a targeted q-profile can be achieved.In order to stay at a sufficiently fast run-time of a few minutes a Bohm/gyro-Bohm model is used for transport.The model does feature free parameters, which have been fitted, utilizing a set of reference discharges.A good fit for T e has been found, while the fit for T i is worse but was found to still be acceptable.Switching to a neural network based version of a higher fidelity model (such as [52]), which would preserve a similar run-time may be interesting once a version applicable to the AT parameter space becomes available.Alternatively, using a more extensive, physics based model with a larger of reference discharges such as [49] might improve reliability.Using any of these options would also be expected to improve the accuracy of the T i prediction.
The model has been used successfully to guide the path to a stable scenario on the 'early heating' approach.While iterative schemes used in between discharges such as ILC [53,54] have been shown to be successful in generating feed-forward trajectories for repetitive control problems, the approach presented here is not employing those.Instead, the model was used in order to generate a setup for the scenario.The result was checked for machine and MHD stability limits and then run, showing good agreement between the predictions and the experimental behaviour.In-between discharge runs in order to do changes to the setup were not done.The model does however have the capability to do so if one wanted to change the target q-profile.Note that while the model is fast enough to operate in between discharges of AUG, high quality q-profile estimates are not available in between discharges.
The resulting scenario can reach a stable bootstrap fraction of close to 50% when operating at β N ∼ 2.7.Above this value, the scenario becomes unstable to NTMs.This threshold is lower than the ideal limit at β N ∼ 3.1 seen in a previous discharge, likely due to a change in transport behaviour after an intermittent drop in heating power and correspondingly T i profiles.Testing whether stable operation can be sustained at higher β, when such an event is not present may be of interest for future experiments.Especially when both off-axis beams are available, achieving fully non-inductive operation in such a scenario seems possible.
A control system based on real-time q-profiles, which is currently under development, may improve NTM stability in the future.In the meantime, optimizing the q-profile in order to avoid low-shear regions around resonant magnetic flux surfaces and thereby improving NTM stability is a promising approach.Using an optimizer for the parameters of interest, such as described in [48] would also eliminate or at least reduce the amount of iterations required to reach a desirable configuration.
To gain information about scenarios in future machines, studying cases with a higher plasma current and thereby lower q 95 would be of interest.Some results for this are presented in [48], further analysis of such a case is planned for the future.
To identify possible additional dependencies in the model, application AUG scenarios with different current and magnetic field is foreseen.Experiments on different machines are planned to test how well the model can be generalized.In the current configuration, the model can in principle be applied to a different scenario or machine if analogous inputs are available.Assuming no required changes, this would mean just changing the input parameters.If inaccuracies were to be found, the behaviour can be adjusted by tuning the free parameters using a new set of reference discharges.However, as there is currently no option to easily do a new parameter fit, considerably more work would be required to do so.Work on this is ongoing and will be presented in a future publication.
Agreement No. 101052200 -EUROfusion).Views and opinions expressed are however those of the author(s) only and do not necessarily reflect those of the European Union or the European Commission.Neither the European Union nor the European Commission can be held responsible for them.

Figure 1 .
Figure 1.Top-down (a) and cross-section (b) plot of AUG showing the paths of the different NBI beams.The magnetic axis and separatrix are shown in blue.

Figure 2 .
Figure2.NBI driven current as calculated by RABBIT[18] in ASTRA[19,20] for each source individually, note that in the real discharge (see section 4.1) not all beams were running.

Figure 4 .
Figure 4. Schematic view of the model setup: ASTRA uses the applied transport model to calculate the temperature and safety factor profiles based on the applied heating systems and underlying plasma density.

Figure 5 .
Figure5.(a) Contributions to the heat conductivity; the blue dashed line shows the behaviour of the gyro-Bohm term towards the edge, a linear approximation towards the edge is used to improve agreement with experimental data.The black dashed line shows the heat conductivity calculated by ASTRA, when using experimental Te-profiles and heating power as input.(b) Effect of the linear approximation in the heat conductivity; the offset in the core between the case with and without the approximation could be mitigated by changing the free parameters; the insert shows a focus at the kink at the transition between edge and core model.

Figure 6 .
Figure 6.Impact of Z eff on q comparing a case at the lower and higher end of the typically expected range.

Figure 7 .
Figure 7. L-H transition threshold compared to power at the separatrix.

Figure 8 .
Figure 8. Conductivity in L-mode.The neoclassical contribution only contributes the nonzero behaviour inwards of ρtor ∼ 0.15.

Figure 10 .
Figure 10.(a) Heating setup in relation to the plasma current, comparing the early heating setup (39221) in the lower plot, with the late heating setup (37722) in the upper plot; (b) distribution of driven current in both scenarios after all heating systems are turned on, the dashed line shows a smoothed version.

Figure 11 .
Figure 11.Time evolution of the safety factor for the late-heating reference scenario, the scenario prediction in ASTRA, and the experimental result of that scenario.The experiment reproduces the behaviour, predicted by the model.

Figure 12 .
Figure 12.Electron temperature time evolution (a) and profiles (b), comparing the experimental results of discharge 39221 with the modelled behaviour.An agreement within the error bars can be seen.

Figure 13 .
Figure 13.Safety factor profiles at some representative time-points, comparing the experimental results of discharge 39221 with the modelled behaviour.A good agreement can be seen.The discrepancy at 2.4 s is caused by different temperature profiles due to fluctuations.

Figure 14 .
Figure 14.Ion temperature time evolution (a) and profiles (b), comparing the experimental results of discharge 39221 with the modelled behaviour.While the general trend is reproduced the general agreement is considerably worse than for the electron temperature.

Figure 15 .
Figure 15.(a) Safety factor time evolution at three different locations comparing experimental data to a simulation with T i fixed to the measured value and fixed to 1.5 times the measured value: (b) simulated β pol compared to the IDE calculation, a very similar behaviour as the core T i can be seen.

Figure 16 .
Figure16.Performance of the validation scenario; starting at the marked time (~7 s) an increase in I nb can be seen without a corresponding change in P NBI due to a change from an on-axis to an off-axis source.

Figure 17 .
Figure 17.Current distribution over time for a discharge, where the NBI power was ramped up until disruption; the scenario becomes fully non-inductive at the end.

Figure 18 .
Figure 18.ECE phase of the mode (a) and mode growth in the magnetic signal (b).

Figure 19 .
Figure19.Current distribution over time for a scenario at higher β N ; At ~6.5 s the discharge switches to a piggyback segment with no relevance to this work, after a locked mode appears.

Figure 20 .
Figure 20.(a) Ion temperature evolution of discharge 41086 at multiple radial locations.While the outer values are mostly unaffected, a significant drop in temperature at t ∼ 5.2 can be seen for the inner values.(b) Logarithmic temperature profiles and Te/T i before (dashed lines) and after (solid lines) the T i drop at ~5.2 s.

Table 1 .
Fitting factors for the free parameters in the model: the values actually used in the model and the mean over the best fitting parameter for each discharge of the reference set; also shown is the range in which the parameters can be altered without significantly changing the quality of the fit; the fit quality of the used values is χ 2 Te = 1.2 and χ 2 Ti = 3.6.