ICT e supercalcolo al servizio di ricerca e imprese

Simulazioni ibride (fluide/Particle-in-cell)
per la fusione termonucleare
G. Vlad, S. Briguglio, G. Fogaccia, V. Fusco
ENEA for EUROfusion, via E. Fermi 45, 00044 Frascati (Roma), Italy
ICT E SUPERCALCOLO
AL SERVIZIO DI RICERCA E IMPRESE
RISULTATI E PROSPETTIVE
17 marzo 2015
ENEA – Via Giulio Romano n. 41, Roma
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
1
Fusion Unit
EUROfusion
Introduction
• 
• 
Controlled thermonuclear fusion is one of the most promising energy sources
for the next near future.
Reproducing in laboratory the nuclear processes which take place in the core of
stars is one of the major challenges of the present day research.
• 
Thermonuclear fusion occurs when light elements (like
Hydrogen or its isotopes) fuse together into new elements,
like Helium, releasing in that process a large amount of
energy.
• 
In order to fuse the light elements together, it is necessary to heat them to
energies of the order of several tens of KeV: in this condition, the gas is highly
ionized. If also high density and good thermal insulation is obtained, the ionized
gas (“plasma”) will undergo a large amount of fusion reactions and the process
will become energetically favourable.
The most promising approach considered by the fusion community is the socalled “Magnetically Confined Fusion”, and the most advanced experimental
devices are the “Tokamaks”.
• 
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
2
Fusion Unit
EUROfusion
ITER
~ 40 m
The next International
Thermonuclear Fusion
experiment (ITER)
The plasma is confined in a
toroidal chamber by a very
high (~6T) magnetic field
ITER is under construction in Saint Paul-lez-Durance
(France). First plasma: ~ 2020; D-T operation: ~ 2027
G. Vlad
Once the plasma, which
curries a strong toroidal
current (~10 MA), is
produced, the topology
of the magnetic field
becomes helicoidal
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
3
Fusion Unit
EUROfusion
Ignition device
• 
• 
• 
• 
• 
• 
The plasma will be heated to the required temperature (~10 KeV) by joule
heating, strongly energetic neutral beams, radio frequency waves, …
Several challenging issues are faced on the way of controlled thermonuclear
fusion:
-  technology (high magnetic fields, superconductor coils, …)
-  materials (heat exhaust, wall loading, blanket, tritium breeding, …)
-  physics (plasma heating, energy confinement, MHD, turbulent transport, …)
Once the thermonuclear reactions become dominant, the plasma temperature will
be sustained by the high energy α particles (Helium nuclei, 3.5 MeV)
It is crucial to have well confined energetic particles (α’s, beams particles, radio
frequency accelerated particles, …) to allow them to slow down and release their
energy by collisions thus heating the bulk plasma
Typical velocity of α particles in an ignited device is of the same order of the
Alfvén velocity (the velocity of propagation of a magnetic field perturbation)
If an electromagnetic perturbation growths in time, because of the resonant
interaction with the energetic particles, the confinement of the energetic particles
themselves can be strongly reduced, before they are able to release their energy to
the bulk plasma, avoiding the “ignition” of the device.
Fusion Unit
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
4
EUROfusion
Theory and modelling
•  The development of theoretical and computational models for the
study of the physical phenomena which determine the plasma
dynamics is of significant importance for the success of future
experiments (“ignition” regime not yet observed in present devices,
extrapolation from present “sub-ignited” regimes is required).
•  We focus our activity in studying the interaction between Alfvén
waves and energetic particles (as, e.g., fusion α’s, beams particles,
radio frequency accelerated particles, …): linear dynamics,
turbulent transport, non linear saturation.
•  The computational model considered is the so-called “hybrid
model”, where a thermal component of the plasma is treated as a
fluid, described by MagnetoHydroDynamics equations (MHD) and
the energetic particles are treated kinetically
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
5
Fusion Unit
EUROfusion
Single particle motion in torus
•  The single particle motion in a tokamak can be very complicated: particles
rapidly gyrate perpendicularly to the equilibrium magnetic field (“gyro-motion”)
while transiting along the torus (“circulating particles”) or experiencing a almost
closed orbits bouncing back and forth (“trapped particles” and “banana orbits”)
because of the characteristic magnetic well of the tokamak configurations;
moreover, those trapped particles experience also a precession motion along the
torus.
•  Thus several characteristic frequencies of energetic particles are present, which
can resonate with the frequencies of the Alfvénic waves, eventually driving them
unstable:
kinetic treatment is important!
circulating particle (“transit frequency”)
Precession of a trapped particle
(“bounce and precession frequencies”)
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
projection of
motion on the
poloidal plane
6
Fusion Unit
EUROfusion
Particle orbit - 1
Trapped particle (“banana” orbit) with precession motion
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
7
Fusion Unit
EUROfusion
Particle orbit - 2
Trapped particle (“banana” orbit) with zero precession
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
8
Fusion Unit
EUROfusion
Particle orbit - 3
Circulating particle
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
9
Fusion Unit
EUROfusion
Alfvén wave (TAE)
frequency spectrum (ω,r):
global mode (TAE)
Alfvén continua
poloidal cross section (R,Z)
G. Vlad
constant flux surface (φ,χ)
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
10
Fusion Unit
EUROfusion
Hybrid model
•  Hybrid model: extend MHD equations by adding a coupling term
with energetic particles (EPs)
• The coupling term between bulk specie and EP
• The EP term needs to keep the details of the
D+
divergenza
tensore pressione
energetiche:
specie
is the divergence
of the EPparticelle
pressure tensor
wave-particles resonances: thus, we need a
kinetic formalism
@%
+ r · (%v) = 0 ,
(9)
• fE is the EP distribution
function described
@t
by Vlasov eq. (collisionless Boltzmann eq.):
dv
1
%
=
rP r · PE + J ⇥ B ,
(10)
dt
c
✓ ◆
d P
@fE
q
= 0,
(11)
+
v
·
r
+
(E + v ⇥ B) · rv fE = 0
r
dt %
@t
m
1
E + v ⇥ B = ηJ,
0,
(12)so-called “gyrokinetic
• solved using the
c
formalism”
1 @B
Z
r⇥E =
,
(13)
c @t
Pi,j;E ⇠ vi vj f (r, v, t) dv
4⇡
r⇥B =
J,
(14)
c
• E(t), B(t) are the solution of the MHD
r · B = 0.
(15) the forces in the Vlasov
eqs. and provide
eq.
d
@
=
+v·r
• vi, vj are the(16)
EP velocity components
dt
@t
Fusion Unit
Workshop
“ICT
e
supercalcolo
al
servizio
di
ricerca
e
imprese
risultati
e
prospettive"
17
marzo
2015
ENEA
Sede
EUROfusion
G. Vlad
11
law:
Computational models
MHD models (field solver)
HMGC:
•  Reduced MHD
•  !two fields;
•  simple geometry;
•  finite difference in r,
•  Fourier in θ !mpol,
•  Fourier in φ ! ntor;
•  fully non-linear.
HMGC
HYMAGYC:
•  fully resistive MHD
•  !14 fields;
•  arbitrary geometry
(curvilinear coordinate
system);
•  finite elements in s,
•  Fourier in χ !mpol,
•  Fourier in φ ! ntor;
•  only linear MHD
dynamics considered.
HYMAGYC
G. Vlad
Gyrokinetic model
•  k ρE << 1 (k is the component of the
wave vector perpendicular to the magnetic
field, ρE is the energetic particle Larmor
radius);
•  magnetic drift orbit widths fully retained;
•  energetic particle pressure: P , P||;
•  solved by particle-in-cell (PIC) techniques;
•  fully non-linear.
• 
• 
• 
• 
• 
k ρE = O(1);
gyro-average fully retained;
energetic particle pressure tensor: Pij;
solved by particle-in-cell (PIC) techniques;
fully non-linear.
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
12
Fusion Unit
EUROfusion
Computational requirements
Goal: ITER relevant case, e.g., n = -30 TAE driven by energetic particles
MHD (field solver) module of HYMAGYC
•  axisymmetric (2D) equilibrium
•  discretization scheme: FE in radius (s), Fourier in
generalized poloidal (χ) and toroidal (φ) angles
•  linear MHD: ! single toroidal mode number: ntor= -30
•  # MHD eqs.:14 eqs.
•  # radial mesh points: ns= 1000;
•  # poloidal Fourier components: mpol= 100;
•  linear system with # eqs.: ns✕mpol✕14 = 1.4✕106
•  matrix elements (double complex): (# eqs.)2 = 1.96✕1012
•  maximum non-zero matrix (block tridiagonal) elements:
456✕ns✕(mpol)2=4.56✕109
•  parallelization of the field solver done in collaboration
with EFDA-HLST using MUMPS (MUltifrontal
Massively Parallel sparse direct Solver), mainly to gain
memory availability
•  runs on CRESCO4 (also on HELIOS (IFERC), Japan)
Gyrokinetic module of HYMAGYC
• 
• 
• 
• 
• 
• 
• 
• 
PIC model;
# particles: np=ns,χ,φ ✕ nppc;
# particles per cell: nppc=512 (v space);
# poloidal mesh points: nχ≈ 8 ✕mpol,max=800;
# toroidal mesh points: nφ≈ 8 ✕ ntor=240;
ns,χ,φ= ns✕nχ✕nφ= 192✕106;
# particles: np= ns,χ,φ ✕ nppc≈ 100 G
! particles memory Mp: 7 double real
variables per particle, Mp=5.5 Tb
•  ! Gyrokinetic module parallelized using
Hierarchical MPI+OpenMP scheme (MPI
inter-node, OpenMP intra-node)
•  typical cpu time/particle/step: 3x10-6 s
•  typical simulation: nsteps=5x104
Memory
cpu time (4576
cores) CRESCO4
(estimate)
cpu time (72000 cores)
HELIOS(estimate)
5.5 Tb
tsimulation≈ 40 gg
tsimulation≈ 2.5 gg
Memory
cpu time (256 cores), MUMPS
cpu time sequential solver
72.96 Gb
tinversion≈ 200 s (inversion)
tbs ≈ 1 s (backsolve per step)
tinversion≈ 3350 s (inversion)
tbs ≈ 29 s (backsolve per step)
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
13
Fusion Unit
EUROfusion
ITER TAE example n = - 30
poloidal cross section (R,Z)
G. Vlad
constant flux surface (φ,χ)
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
14
Fusion Unit
EUROfusion
JET model equilibrium
• 
• 
• 
• 
n=4 non-linear simulation,
soliton avalanches,
energetic particles displaced toward the
outer edge of the torus,
critical phenomenology,
! possibly preventing ignition!
n4_soliton_avalanches.m4v
movie caption:
Electrostatic potential
ϕ Fourier
components vs r
energetic particle
pressure radial profile
Electrostatic potential
ϕ structure in (R,Z)
frequency spectrum
(ω,r):
Energetic particle
driven mode (EPM) +
Alfvén continua
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
15
Fusion Unit
EUROfusion
Non-linear dynamics studies
Energetic electron driven internal kink instability: “e-fishbone”:
•  MHD mode driven by energetic electrons (observed in present devices, similar dynamics of
certain energetic-ion instabilities driven in ignited plasmas)
•  HMGC suited for detailed studies of non-linear saturation by Hamiltonian mapping techniques
•  plots of energetic particles in the plane (Θ,Pϕ), with Θ the wave phase seen by the energetic
particles and Pϕ the toroidal angular momentum.
peaked-off_eps0.1_pphi_phase_den_power_res_grande-desktop.m4v
G. Vlad
Workshop “ICT e supercalcolo al servizio di ricerca e imprese risultati e prospettive" 17 marzo 2015 ENEA - Sede
16
Fusion Unit
EUROfusion