NIMROD Team Meeting Minutes

May 1-2, 2009

Boulder, Colorado

Performance Development and Scaling

Carl Sovinec described the implementation of Newton-Krylov iteration in NIMROD’s two-fluid implicit leapfrog advance. This iteration centers the algebraically nonlinear momentum-density advection term in the velocity advance and the nonlinear Hall term in the magnetic field advance. ‘Algebraically nonlinear’ here means nonlinear with respect to the change over each time-step; most of the physically nonlinear contributions are explicit. The temporal staggering of the implicit leapfrog and the bilinear character of the two terms facilitate computation of the Jacobian for each Newton step, and the implementation is a modification of previously existing code. The nonlinear residual also has a simple form with the bilinear nature of these terms, and iteration tolerance is checked with the norm routine from NIMROD’s 3D linear solver. Initial tests of the nonlinear iteration do not show significant changes, but the algorithm may improve nonlinear stability.

Jeong-Young Ji gave a presentation on closure theory and emphasized that use of a transport equation to close a moment equation is not consistent. In the general moment approach, all equations have the same general form. Terms associated with derivatives of the background Maxwellian distribution lead to drives in a matrix equation for the moments of the kinetic distortion. He has considered Braginskii (small gyro-orbit), neoclassical (high to low collisionality with flux-surface averages over assumed nested surfaces), and general closure without flux-surface averaging and with the full linearized Coulomb collision operator. Recent work includes a simple fit for the kernel function that agrees with Chang and Callen in the collisionless limit. The simple fit appears to be very accurate and should help make nonlocal closure computations more practical.

Travis Austin described algebraic multigrid (AMG) and its possible application to NIMROD. The motivation is scalability as MG has been the only matrix algorithm to show O(N¹) performance for Poisson problems on the largest parallel machines available today. The geometric version determines restriction and prolongation operations geometrically, but this has problems with stretched meshes. Algebraic versions determine strongly connected points in the matrix graph and use the strong connections to define coarse approximations. This is challenging for high-order finite elements, because the matrices are not M-matrices. His approach is to use linear elements to precondition, based on ideas from Paul Fischer. Important algebraic modes within each high-order element are determined by separate minimization problems. This form of problem-reduction could also be used with direct solvers when they are applied to preconditioning.

Jianhua Cheng presented a hybrid algorithm that uses delta-f on the full kinetic equation for ions. It adds the delta-B-parallel response and is expected to be useful for edge physics and for the MHD ordering. The electrons presently use fluid modeling, assuming an isothermal response for now, and displacement current is neglected. The algorithm is implicit for compressional-wave physics at high beta, and Cheng has made the advance second-order accurate using iteration for the particle/electric field coupling. Test problems considered so far include compressional and shear Alfven waves, the ion acoustic wave, and the Whistler wave. Work on Harris-sheet reconnection is starting. Modifications for drift and gyrokinetic electrons are also being considered.

Project Improvement

Scott Kruger gave a presentation and a demonstration of regression testing, which uses the txtest program written by David Wade-Stein. Its main purpose is to ensure that development in one area does not affect other parts of the code. It is used for both automated testing, which is run nightly, and manually by developers as they are working. Test cases should be short and coarsely resolved, as tests that take a long time discourage use. However, this means that we may be looking for changes at high precision (to a controllable tolerance via h5diff), so deciding when a discrepancy is acceptable or not can be tricky, as in the ORNL tearing-mode case. More tests are needed, including cases from Schnack’s V&V list. The test cases are placed in separate SVN repositories that are organized by machine and require 1) a shell script to run the case, 2) input files, and 3) an edit to makefile.am. There is a nimrod-regress mailing list for those who would like to see all output.

Applications

Charlson Kim described simulation and visualization work for ICCs at the PSI-Center. He has been simulating the Caltech spheromak/plasma jet experiment. The simulated injector specifies toroidal field as a boundary condition and has an insulating layer above it. Visualization uses the Visit package from LLNL. It is used to show the formation of closed-flux regions during decay. It has also been used to show 3D visualization capability.

Ryoji Takahashi gave a presentation on kinetic effects of fast particles with resistive MHD instability. Results from JET suggest that large hot-particle beta breaks the empirical rho-star scaling of the 2/1 mode from smaller machines. Linear computations show that the stabilizing effect of hot particles is strongest in the resistive regime (with respect to the ratio of beta-normal and internal inductance) of this mode, so it’s not an ideal resonance effect. The study will first try to explain the linear effects then move to nonlinear. However, the lack of nlayer parallelization with particles is a hindrance. An action item is for Takahashi and Kim to first review the implementation to verify that characteristics and weight-evolution equations are correct for nonlinear conditions.

Dalton Schnack discussed efforts to model giant sawteeth in DIII-D. The work involves collaboration with Alan Turnbull at GA and seeks to model the effects reported in the PoP paper by Choi. Fast particles stabilize the internal kink, but slow changes in the magnetic shear lead to less stabilization over the cycle. Thus, the crash results from a changing kinetic effect and not from changes to MHD stability; though, the hot-particle beta also changes. The high-energy tail is actually from RF kicking beam-driven particles to 100s of keV, and these particles are the important ones according to Turnbull. The NIMROD hot-particle equilibrium presently does not model this part of the distribution function, and extremely energetic particles will lead to very small timestep. [Brute forcing with better parallelization may be possible, at least linearly.] Schnack noted that the analytical relation for kink stabilization/fishbone destabilization seems to require two-fluid effects. Also, NIMROD’s two-fluid advance appears to be bypassed completely when particles are used. A list of action items for this study is:

1. Confirm that we are doing an appropriate and tractable case. (Schnack)

2. Correct particle acceleration and field advance for running particle with 2-fluid. (Brennan, Takahashi, and Kim)

3. Modify the distribution function for an RF tail. (Kim)

4. Improve particle parallelization, possibly through Eric Held’s code for closures. (Kim, Kruger, and Held)

Tom Jenkins described his progress on 2/1 island saturation subject to a localized (RF-like) source term. If it is placed in the right location in the poloidal plane, the localized source can completely stabilize the island. However, rational surfaces move in response to the source. This is an important effect that warrants a more thorough investigation than it was given in the Pletzer-Perkins paper. Another consideration is the prompt response of the profile (after several Alfven times) vs. the long-time response after resistive diffusion times when the source drives net current. The former is likely more important for active steering systems that are being used in experiment. Jenkins also discussed progress in coupling with the GENRAY ray-tracing code through use of the ‘plasma state’ that has been defined as part of the SWIM project.

Valerie Izzo described her progress in simulating runaway electrons during tokamak disruption. In her modeling, the seed electrons are nonrelativistic, and disruption is induced by a modeled gas jet. The electron trajectories are now computed at runtime through the use of code developed by Eric Held for nonlocal closure, instead of separately as postprocessing. The computations compare profiles with low and high elongation; the high elongation cases have more nonlinear coupling. Next steps will include perpendicular drifting, so that orbits need not be perfectly tied to field-lines, and synchrotron radiation effects.

Ping Zhu gave an update on his two-fluid ballooning computations. With some unexplained fix of memory use on Franklin, he is again able to run linear computations with high resolution. From equilibria generated from ESC, NIMROD produces ideal ballooning modes with little dissipation across most of the spectrum. However, Zhu finds that particle hyperdiffusivity alone does not provide enough dissipation to avoid numerical instability for all n; some thermal conduction is also important. At very low n, the computations can still produce slowly growing noise in n in the transition region between small and large elements. With two-fluid effects in Ohm’s law and gyroviscosity, the modes show up-down asymmetry. In these cases, Zhu is not expecting FLR stabilization, just a reduction in growth rate. Here, the very low n behavior seems to be okay.

Bonita Burke (formerly Squires) reviewed verification work for ideal ballooning behavior through a benchmark with GATO and ELITE. Her recent sets of computations have no resistivity gradient in the pressure pedestal, and mass density has a hyperbolic-tangent profile. Burke emphasized the need to scan S-values in the core and halo regions separately for each case. The existence of a halo is most important for n-values of around 5 or 6. For the equilibria labeled dens8 (the more unstable case), NIMROD shows good quantitative agreement with GATO for the n-values where it has been applied (n>11). ELITE’s predictions are close and depend on whether parallel compression energy is used. Burke reported that NIMROD also finds the correct point of marginality between n=8 and n=9 for the less unstable dens6 case (GA did not provide their spectra); though, ‘false’ modes with small non-convergent growth-rate can appear in computations below marginality. A new equilibrium is being used for peeling/ballooning tests, and Burke has computed an ideal spectrum. However, with a realistic Spitzer resistivity profile, a peeling peak is not evident, possibly due to faster resistive ballooning at high n.

Mark Schlutt gave an update on his application of NIMROD to straight stellarators. The objective of the study is to determine what limits beta, how it occurs, and what influences gradual magnetic topology change. An important aspect is the pressure-induced Pfirsch-Schluter current. The NIMROD results will be compared with analytical work for 3D islands with finite parallel thermal conductivity. Schlutt is considering two different problems, one that applies nonuniform heating to a 2D magnetic configuration and the other that adds 3D magnetic perturbations. The computations start with nested flux surfaces, and induced resonant harmonics change the topology. Schlutt compared low and high heating power, and high heating produces poloidally elongated configurations that can split into m=2, n=2 islands.

Sovinec gave an overview of helicity-injection related computations by students Tom Bird, John O’Bryan, and Eric Howell. Bird had presented his work at the 2009 U-WI Undergrad Symposium. His simulations of current injection and relaxation in the Pegasus ST use the toroidally localized version of the ad hoc source used in Jenkins’ study. The latest shows helical current paths distorting from self-induced field overcoming vacuum field, but full merging and current multiplication have yet to occur with a basic resistive-MHD model. New computations are using temperature-dependent coefficients. O’Bryan is presently performing 2D modeling of flux compression from an assumed, relaxed state. The computations use T-dependent coefficients with magnetization and Ohmic heating. Computations with realistic parameters are showing more heating than current drive. Howell is studying the linear properties of decaying spheromak profiles, which can now be generated systematically with NIMEQ. He is first checking the S-scaling of resonant modes with resistive MHD and will then consider the influence of two-fluid effects.

Jacob King reported on progress with two-fluid tearing and dynamo. The present set of computations provide scans over normalized sound gyroradius with uniform pressure to avoid equilibrium diamagnetic effects for now. Single-helicity cases started from paramagnetic pinch profiles produce a saturated island width that is approximately independent of sound gyroradius. The out-of-phase components of the magnetic perturbations do change, however, and this is important for the Hall dynamo contribution. The dynamo spatial structure is more complicated in the two-fluid cases and shows significant cancellation of the Hall and MHD dynamo contributions away from the rational surface. The two-fluid cases also show slow rotation after saturation. King also described some preliminary results for single-helicity axis states produced in RFX. The modeling will investigate the importance of temperature-dependent transport effects.