Physics Research in The SciDAC Center for Wave – Plasma Interactions

1. Physics Research in The SciDAC Center for Wave – Plasma Interactions Paul Bonoli on behalf of the SciDAC Center for Simulation of Wave-Plasma Interactions 17th Topical Conference on Radio Frequency Power in Plasmas Clearwater, Florida May 7-9, 2007

2. Participants in the Center for Simulation of Wave – Plasma Interactions L.A. Berry, D.B. Batchelor, E.F. Jaeger, E. D`Azevedo, M. Carter D. Smithe C.K. Phillips, E. Valeo N. Gorelenkov, H. Qin P.T. Bonoli, J.C. Wright R.W. Harvey, A.P. Smirnov N.M. Ershov M. Brambilla R. Bilato Politecnico di Torino R. Maggiora V. Lancellotti M. Choi D. D’Ippolito, J. Myra - Lodestar Research

3. What is “SciDAC” ? • Scientific Discovery Through Advanced Computing: • An initiative of the US DoE to take advantage of high performance computing platforms (Massively parallel systems) to conduct physics research. • This initiative brings applied mathematics and numerical algorithms expertise to the physicists. • Our project also has a significant theoretical and experimental component that has been critical to its success.

4. The role of RF power in fusion plasmas has evolved over the years • Early applications involved bulk plasma heating and non-inductive maintenance of the entire plasma current: • Ion cyclotron resonance heating (ICRH), lower hybrid heating (LHH), and electron cyclotron resonance heating (ECRH). • Lower hybrid current drive (LHCD) was especially successful and efficient. • More recent applications have utilized RF waves for localized control of the plasma current and pressure profiles: • LH and EC current drive for sawtooth control, NTM control. • Mode converted ICRF waves for shear flow generation and current generation. • Current profile control for access to high confinement (‘advanced tokamak”) regimes with high bootstrap current fraction (> 70%). • Most of the RF applications listed above are planned for ITER: • A predictive capability is needed to insure success in ITER applications.

5. So what are some of the challenges ? • Complete description of the wave-particle interaction involves integrating several “separate” calculations: • Antenna coupling (linear and nonlinear response) • Wave propagation • Wave absorption • Multiple spatial scales can occur: • ICRF Mode conversion and current drive • RF sheath formation in the edge • 3D reconstruction of RF fields and RF power deposition is required when considering coupling of RF response to MHD and turbulence models. • RF waves can interact in a nonlinear fashion with thermal and nonthermal particles in the plasma: • Nonthermal electron and ion distributions produced by the RF wave itself. • Nonthermal ions already present in the plasma due to neutral beam injection (NBI) and fusion reactions (alpha particles).

6. Two full-wave solvers have been advanced within the RF SciDac Initiative Blowup region Slow ion cyclotron wave • All Orders Spectral Algorithm (AORSA) – 1D, 2D & 3D (Jaeger) • Spectral in all 3 dimensions • Cartesian/toroidal coordinates • Includes all cyclotron harmonics • No approximation of small particle gyro radius r compared to wavelength l • Produces huge, dense, non-symmetric, indefinite, complex matrices • TORIC – 2D (Brambilla/Bonoli/Wright) • Mixed spectral (toroidal, poloidal), finite element (radial) • Flux coordinates • Up 2nd cyclotron harmonic • Expanded to 2nd order in r/l • Sparse banded matrices Electrostatic ion Bernstein wave

7. Wave propagation and the plasma response are governed by the Maxwell-Boltzmann system of equations For time harmonic (rapidly oscillating) wave fields E with frequency ω, Maxwell’s equations reduce to the Helmholtz wave equation: Wave Solvers (AORSA) (TORIC) The plasma current (Jp) is a non-local, integral operator (and non-linear) on the rf electric field and conductivity kernel: The long time scale response of the plasma distribution function is obtained from the bounce averaged Fokker-Planck equation: Plasma Response (CQL3D) 0 where Need to solve this nonlinear, integral set of equations for wave fields and velocity distribution function self-consistently. This requires an iterative process to attain self-consistency.

8. Calculation for C-Mod minority H, NR = 128, NZ = 128,[256 processors for 3 hrs on Cray XT3 – ORNL]

9. H+ ~0.1% B Experimental measurements of the energetic ion tail on C-Mod have been made using a compact neutral particle analyzer (CNPA) – V. Tang • Using a Maxwellian fit to data gives Tion  70 keV. • Vincent Tang, MIT [PhD Thesis, 2006; also V. Tang, PPCF, 2007] • Good agreement between simulation and measurement on Tion

10. 3D (r, V , V//) distribution function from CQL3D – AORSA reproduces CNPA measurements using a synthetic code diagnostic Building the synthetic diagnostic is a lot of work ! V. Tang, PhDThesis, MIT (2006). Also PPCF, (2007).

11. A Puzzle arose from initial simulations of high harmonic fast wave (HHFW) – fast ion beam interaction in DIII-D Stronger Beam Interactions at 4D (60 MHz) Than at 8D (116 MHz) Observed in DIII-D CQL3D-AORSA predicted increased absorption as frequency was raised. Monte Carlo ORBIT code (ORBIT-RF) combined with an RF operator (using fields from TORIC solver) did reproduce the experimental trend. Could finite ion orbit effects be playing a role ? • DIII–D high density L-mode 1014 (/s) Neutron reaction rate Power Sn: neutron enhancement factor 217-05/MC/jy

12. AORSA and CQL3D were iterated to solve for the wave fields and distribution function self-consistently in the DIII-D HHFW – NBI Interaction f = 60 MHz (4th harmonic D; non-Maxwellian) and 2% minority H (2nd harmonic H, Maxwellian) 0th iteration PRF = 0 1st iteration PRF = 1.1 MW 5th iteration PRF = 1.1 MW 6th iteration PRF = 1.1 MW CQL3D: Neutral beam D 2Ω(H) Prf 4Ω(D) 4Ω(D) 4Ω(D) 4Ω(D) 2Ω(H) e 2Ω(H) 2Ω(H) e e e H P(D) = 0.29 MW = 23.8 % P(D) = 0.84MW = 66.3 % P(D) = 0.76 MW = 58.2 % P(D) = 0.75 MW = 57.1% At 60 MHz, about 57% of the RF power is absorbed by the deuterium

13. Iterative solution for non-Maxwellian H and D in DIII-D at 60 MHz (2% H assumed) shows significant absorption of HHFW on background H

14. We are investiagting finite ion orbit effects using two approaches: • The diffusion coefficient (D) has been evaluated by a direct orbit integration using electric fields from AORSA: • The “DC” code calculates the bounce average RF (not quasilinear) diffusion coefficients as a function of (v//, v, r) • The Monte Carlo code ORBIT RF has been combined with the TORIC ICRF solver: • Self-consistent iteration not yet carried out.

15. Orbit widths not negligible, even in C-Mod • Shown at right are trajectories for 12 particles in the C-Mod case: • All particles launched from 80 cm • 4 equi-spaced || velocities • 3 equi-spaced  velocities • 409,600 complete poloidal orbits • Diffusion Coefficient calculations done on CRAY XT3 (ORNL) using 256 processors @ 10 min.

16. Comparison of RF Diffusion Coefficients from AORSA/CQL3D and Lorentz Force Integration at R= 72 cm(C-Mod minority ICRH case) • Lorentz Orbit Code -DC • AORSA/CQL3D

17. ORBIT-RF has been used to show that large RF electric field can destroy superadiabaticity via phase stochasticityMechanism works even in the absence of collisions Using parameters from the C-Mod minority heating case: PRF = 0.6 MW ERF 1 kV/m -> =0.428 • Mapping parameter () from Stix is:  Y/a Green = 0.1 kV/m Red = 1 kV/m X/a  

18. Impact on experiments – TORIC and AORSA are now being used to design RF flow drive experiments TORIC E field is used in synthetic PCI on C-Mod that measures the perturbed density due to mode converted ICRF waves. Y. Lin, A. Parisot, J. Wright; PoP, 2005, PPCF, 2005 • Shear flow generation is possible in fast wave to ion cyclotron wave mode conversion (Myra, Jaeger et al Phys Rev Lett, 2003 and Melby et al Phys Rev Lett. 2003)

19. Comparison of TORIC (FLR) and AORSA (all orders) fieldsMode Converted Wave Fields are Similar ICW TORIC at 240Nr x 255 NmAORSA at 230Nx x 230 Ny • Both codes are using the same equilibrium from an Alcator C-Mod discharge with mixture of D-3He-H in (21%-23%-33%) of ne proportion. IBW FW

20. Calculations on the Cray XT-3 have allowed the first simulations of mode conversion in ITER ITER with D:T:HE3 = 20:20:30 with NR = NZ = 350, f = 53 MHz, n =2.5x1019 m-3 (4096 processors for 1.5 hours on the Cray XT-3) E_perp Blowup (E_parallel) Mode converted Ion Cyclotron Wave (ICW)

21. Capability to efficiently compute 3D wave fields will be important for assessing antenna – edge interaction and for coupling to MHD simulations NSTX ITER Scenario 2 ITER Simulation using AORSA (Fields Summed over 169 toroidal modes) (Jaguar: 2048 processors for 8 hrs)

22. Lower Hybrid Current Drive (LHCD) Simulations in Present Day Devices and ITER • Combined Fokker Planck – ray tracing models are already run within time dependent simulations (DELPHINE-CRONOS, TSC – LSC, TRANSP – LSC, TASK) • Want to incorporate LH full-wave description with Fokker Planck codes (diffraction, wave focusing). • Full-wave simulations can be done now in mid-size present day devices (even including the full launcher spectrum). • Full-wave LH simulation for ITER is a petascale problem.

23. Comparison Between Experimental Hard X-rays and Synthetic Diagnostic (CQL3D) for Shot 1060728104 • HXR Spectra from Synthetic Code are narrower than measured profiles. • This suggests that radial diffusion effects on fast electrons or diffraction could be important. • Next step is to include a model radial diffusion operator in the CQL3D simulations and couple CQL3D to a full-wave solver.

24. Full-wave solver TORIC is now being coupled to the CQL3D Fokker Planck solver • Electron and ion plasma response in TORIC can now beevaluated using the nonthermal fe, from CQL3D(E. Valeo, C.K. Phillips, J. Wright). • Final step is to reconstruct DQL from TORIC and pass it to CQL3D. • Full-wave LH Simulation first performed on the MARSHALL cluster at MIT(6 days @ 32 processors) revealed significant spectral broadening due to diffraction. • Now simulation can be done on CRAY XT3/XT4 JAGUAR at ORNL (1 hour @ 4096 processors). • Coupled CQL3D – TORIC simulation with iteration requires about 33,000 CPU hours per toroidal mode.

25. Alcator C: Full-wave LH Simulation Using TORIC-LH with Maxwellian electron response Field pattern shows evidence of weak, multi-pass absorp- tion

26. Alcator C: Full-wave LH Simulation Using TORIC-LH with quasilinear electron response Field pattern shows evidence of strong absorption on first pass

27. Antenna Coupling Code (RANT, TOPICA) sheath BC Edge Power Dissipation Wave Propagation Code (AORSA, TORIC) Parametric Decay RF Power Coupling to Core Long Range Goal: Integration of linear and nonlinear rf physics in the edge plasma

28. Linear Coupling of an ICRF Antenna with the Edge Plasma can be Computed by Integrating the TOPICA and TORIC Codes • TOPICA: • Fully 3D solid antenna structure model (including FS, box,…) • Parallel code has been used to model the ITER ICRF antenna). Alcator C-Mod: “E” – antenna

29. TOPICA is coupled to TORIC through a surface admittance matrix • Shown is the perpendicular response to a perpendicular drive electric field at the surface.

30. Two Approaches are Being Pursued to Study the Nonlinear ICRF antenna – edge Interaction • Implementation of RF sheath boundary conditions in full-wave solver: • Modify metal wall BC in TORIC solver • Approach will be “spectral” • Time domain simulations using 3D EM field solver - VORPAL • Fully implicit time domain dielectric response module has been implemented for electrons and ions. • Can also use PIC-treatment for ion response (fully nonlinear).

31. Dielectric for the RF Time Domain model has been tested on a cold plasma ICW mode conversion case ICW Full-wave Solution (TORIC) Time Domain Solution IBW

32. VORPAL Time Domain Simulation of Antenna Loop

33. VORPAL Time Domain Simulation of Antenna Loop Surface E-field on Loop antenna RF B-field Surface E-field Wavefronts Radiation Pattern Radiation from Behind

34. TORIC and AORSA can address physics effects that are missing in other fast particle mode simulation tools • Most existing simulation tools for fast particle modes are based on linear or nonlinear MHD / hybrid models that neglect or treat approximately: • FLR effects, finite c effects, cyclotron resonances • Full toroidal geometry, fast ion distributions • TORIC and AORSA retain these effects and, in addition, include much higher spatial resolution of modes, but in a linear or quasilinear treatment • Approach >>> Utilize TORIC to study dynamics of driven modes in the linear regime • begin in the c regime • modify code as needed for c regime • validate code with results from driven mode experiments

35. Compressional Alfven Eigenmodes (CAE’s) were first predicted by the HYM Code (E. Belova) and then Observed in NSTX TORIC Simulation HYM RD IBW CAE? ICW? Nm=511 Nr=483 BT = 0.426T Ip= 0.6 MA qa=12 f = 2.5 MHz

36. Summary • Coupled full-wave Fokker Planck solvers (CQL3D-SIGMAD-AORSA) have been used to simulate minority ICRF, fast wave current drive, and mode conversion current drive in present day tokamaks and in ITER: • Full coupling uses self-consistent nonthermal ion distributions and quasilinear diffusion coefficient in differential form. • Comparison of synthetic diagnostics (CNPA)with experiment is encouraging. • Both AORSA and TORIC can simulate ICW/IBW mode conversion in present day tokamaks and in ITER: • Comparison of synthetic diagnostics (PCI) with experiment is encouraging. • Flow drive has been computed for C-Mod and ITER using AORSA. • Full-wave LH field simulations have been performed for the first time ever in a mid-size tokamak: • Self-consistent coupling to an electron Fokker Planck solver is now underway.

37. Summary • Fully self-consistent coupling of our full-wave solvers to a Monte Carlo orbit code (ORBIT RF) is underway: • Full coupling will use statistical nonthermal ion distributions and quasilinear diffusion coefficient in differential form • We have used the ORBIT RF code to verify the applicability of quasilinear theory in present day ICRF tokamak experiments. • A self-consistent treatment of the RF antenna – edge plasma is underway: • Linear antenna coupling problem is substantially completed using the TORIC – TOPICA suite. • Boundary conditions for sheaths will be implemented in our full-wave solver (TORIC) • Proof of principle time domain simulations of sheath formation have been done using VORPAL. • We have investigated novel applications of our full-wave solvers to simulate driven modes at ci in NTSX.