A generalized weighted residual method for RFP plasma simulation

1. A generalized weighted residual method for RFP plasma simulation Jan Scheffel Fusion Plasma Physics Alfvén Laboratory, KTH Stockholm, Sweden

2. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

3. Basic idea Time differencing numerical initial value schemes (even implicit) require extremelymany time steps for problems of physical interest, where there are several separated time scales. Causality is already embedded in the governing PDE’s – - there is no need to mimic causality by time stepping. Spectral methods (solution expanded in basis functions) are successful in the spatial domain – why not employ them also in the time domain? By expanding in time + physical space + physical parameters, the computational result will be semi-analytical. (Analytic basis functions with numerical coefficients). Ideal for scaling studies, for example.

4. What is the Generalized Weighted Spectral Method (GWRM) ?  Fully spectral weighted residual method for semi-analytical solution of initial value partial differential equations. • All time, spatial and physical parameter domains are represented by Chebyshev series, enabling closed and approximate analytical solutions. The method generalises earlier spatially spectral, finite time difference methods.  The method is acausal and thus avoids time step limitations.  The spectral coefficients are determined by iterative solution of a nonlinear system of algebraic eqs, for which a globally convergent semi-implicit root solver (SIR) has been developed. • Accuracy is controlled by the number of included Chebyshev modes. • Efficiency is controlled also by the use of temporal and spatial subdomains. • Intended for efficient solution of nonlinear initial value problems in fluid mechanics and magnetohydrodynamics, including simulation of multi-time-scale RFP confinement and transport.

5. The Generalized Weighted Spectral Method (GWRM) Consider a system of parabolic or hyperbolic initial-value PDE’s, symbolically written as D is a nonlinear matrix operator, f is a forcing term. D and f contains both physical variables and physical free parameters (denoted p). Initial u(0,x;p) + (Dirichlet, Neumann or Robin) boundary u(t,xB;p) conditions. Integrate in time: Solution u(t,x;p) is approximated by finite, first kind Chebyshev polynomial series. Definition: Chebyshev polynomial Tn(x) = cos(n arccosx). For simplicity – here single equation, one spatial dimension x, one physical parameter p.

6. The Generalized Weighted Spectral Method (GWRM) The Weighted Residual of the GWRM is given by with The TP-WRM coefficients are now obtained from the nonlinear system of algebraic equations where The initial state is expanded as

7. The Generalized Weighted Spectral Method (GWRM) • COMMENTS • • Boundary conditions are transformed into Chebyshev space (using Chebyshev interpolation); • they enter at the highest modal numbers of the spatial Chebyshev coefficients. • • All computations are in Chebyshev space. • • Efficient procedures for integration, differentiation and nonlinear products in Chebyshev • space have been developed. • • Chebyshev polynomial expansions have several desirable qualities: • - converge rapidly to the approximated function • are real and can be converted to ordinary polynomials and vice versa • minimax property - they are the most economical polynomial representation • can be used for non-periodic boundary conditions

8. Simple GWRM example - the linear diffusion equation Solution to be determined: The coefficients aqrs are determined by iterations, using a root solver. for 1≤ q ≤ K + 1 with Boundary conditions enter here

9. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

10. ODE example Light a match – a model of flame propagation: y – flame radius green - exact solution d = 0.05; # Chebyshev modes K = 20, # time domains Nt = 1, error = 0.01

11. ODE example, cont’d • = 0.01 • # Chebyshev modes K = 8 • # time domains Nt = 10 • error = 0.01 • = 0.0001 – Stiff problem! • # Chebyshev modes K = 5 • # time domains Nt = 100 • error = 0.1 • Adaptive grid should be used • for improved accuracy.

12. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

13. SIR - a globally convergent root solver The GWRM - well adapted for iterative methods for two reasons: 1) Basic Chebyshev coefficient equations are of the standard iterative form 2) Initial estimate of solution vector can be chosen sufficiently close to the solution by reducing the solution time interval Instead of using direct iteration, the Semi-Implicit Root solver (SIR) finds the roots to the equations or, in matrix form

14. SIR - a globally convergent root solver The system has the same solutions as the original system, but contains free parameters in the form of the components of the matrix A. The parameters can be chosen to control the gradients of the hypersurfaces . Adjusting these parameters, global, quasi-monotonous and superlinear convergence is attained. In SIR, whereas are finite and is controlled; it produces limited step lengths, quasi-monotonous convergence; and approaches zero after some initial iterations. Newton’s method is a special case of the present method, when all • Rapid second order convergence is generally approached after some iteration steps. • Relationship to Newton’s method - approximately similar numerical work; inversion of a Jacobian matrix at each iteration step.

15. Newton’s method 2D example solution

16. Newton’s method with linesearch 2D example local minimum f ≠ 0

17. SIR 2D example finds solution

18. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

19. Accuracy - the Burger equation • Burger’s nonlinear equation (parabolic) Parameters: Solution compared to Lax-Wendroff (explicit time differencing): GWRM parameters - (S,M,N) = (2,7,6), 13 iterations. L-W marginally stable parameters - 10-3 TP-WRM solution, using two spatial subdomains TP-WRM solution error, as compared to exact solution Result: GWRM 50 % faster than L-W for same accuracy.

20. GWRM Burger equation solution, including viscosity dependence u = u(t,x;v) Initial conditionj(x) = x(1 - x) and boundary condition u(t,0;v) = u(t,1;v) = 0. Solution shown versus x and v at time t = 2.5. Here K = 8, L = 10, and M = 2.

21. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

22. Efficiency Wave equation, forced (hyperbolic) Parameters: Exact solution GWRM solution (averages out fast time scale) (slow + rapid time scale)

23. Efficiency • Forced wave equation solutions u(t,x0) for fixed x = x0 GWRM (K,L) = (6,8) Lax-Wendroff, explicit ∆x = 1/30, 900 time steps Crank-Nicholson, implicit ∆x = 1/30, 100 time steps

24. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

25. Discussion GWRM work so far: • The time- and parameter-generalized weighted residual method, J. Scheffel, 2008. (GWRM method outlined) • Semi-analytical solution of initial-value problems, D. Lundin, 2006. (Resistive MHD stability of RFP and z-pinch) • Application of the time- and parameter generalized weighted reidual method to systems of nonlinear equations, D. Jackson, 2007. (Navier-Stokes equations, Rayleigh-Taylor instability) • Further development and implementation of the GWRM A. Mirza, ongoing Ph D studies (Application of GWRM to nonlinear resistive MHD) SIR: • Solution of systems of nonlinear equations, a semi-implicit approach, J. Scheffel, 2006. (SIR outlined) • Studies of a semi-implicit root solver, C. Håkansson, M Sc Thesis. (Efficient SIR compared to other methods)

26. Discussion (cont’d) • The GWRM is shown to be accurate for spatially smooth solutions - convergence including sharp gradients should be further studied. • Efficiency is central; SIR involves Jacobian matrix inversion of Chebyshev coefficient eqs - N eqs takes O[N3] operations. • Methods to improve efficiency have been developed - temporal subdomains and spatial subdomain techniques using overlapping domains. • Further benchmarking of efficiency for MHD relevant test problems should be carried out as well as corresponding comparisons with implicit methods. • SIR efficiency and linking to GWRM is presently being optimized.

27. OUTLINE • What is the GWRM? • ODE example • SIR - a globally convergent root solver • Accuracy • Efficiency • Discussion • Conclusion and prospects

28. Conclusion and prospects • A fully spectral method, the generalized weighted residual method (GWRM), for solution of initial value partial differential equations, has been outlined. • By representing all time, spatial and physical parameter domains by Chebyshev series, semi-analytical solutions can be obtained as ordinary polynomials. (“Semi-analytical”: expansion in basis functions with numerical coefficients.)  Computed solutions thus contain r- t- and parameter dependence explicitly.  The method is global and avoids time step limitations. • Spectral coefficients are found by iterative solution of a linear or nonlinear system of algebraic equations, for which an efficient semi-implicit root solver (SIR) has been developed. • Accuracy is explicitly controlled by the number of modes and subdomains used.  To improve efficiency, a spatial subdomain approach has been developed.  Problems in fluid mechanics and MHD will be addressed. • Future applications involve studies of nonlinear plasma instabilities at finite plasma pressure in stochastic magnetic field geometries, in particular operational limits in reversed-field pinches.

29. Thank you for your attention!