## Simulating CR particle acceleration in shocks modified by CR nonlinear effects

Jones and Kang (2003) developed a new, fast numerical scheme for the CR diffusion convection equation to be applied to the study of the nonlinear time evolution of CR modified shocks for arbitrary spatial diffusion properties. To reduce the effort required to solve the diffusion convection equation during diffusive shock acceleration, they take advantage of the expected smooth structure of the particle distribution function, f(p), and used a finite volume approach in momentum space along with a simple model for the distribution inside individual volumes of momentum space. The method extends that in Jones (1999), and Jun and Jones (1999), and is similar to the scheme introduced by Miniati (2001), but computationally simpler. Those authors actually ignored spatial diffusion, so could not explicitly treat diffusive shock acceleration. They applied analytic test-particle solutions of the diffusion convection equation forfp) at shock jumps.

The method used in Jones and Kang (2003) is as follows. Ignoring momentum diffusion they write in 1D-geometry:

dt dx dx ^ dx J 3 ^ dx J dp where S is a convenient source term, and all the other symbols take their usual meanings. The number of CR particles in the momentum bin Ap;- = [,p.+1] is pi+1 2 . , n = J p2f(p)p. (4.23.8)

Integrating over the finite momentum volume bounded by Api gives

^ + = Fn - Fn + A\Kn ^LI + S„., (4.23.9) dt dx n' n'+1 dx { n' dx J n'

where

• , , c ,, - p du with the momentum speed , p =---,

3 dx

Pi+1 df 2 lp'+1 df 2 pi+1 df 2 / Kni = J Kfp2dp J fp2dp ^ J Kfp1dp ni, (4.23.11)

Pi+1

The first pair of terms on the right hand side of Eq. 4.23.9 represents the net flux of CR across the boundaries of the individual momentum bins owed to adiabatic flow compression or expansion. Extension of this term to include fluxes owed to other energy loss or gain processes, such as momentum diffusion or radiative losses, is obvious. In addition there was defined

Pi+1

For relativistic particles g, is proportional to the CR energy in the bin. With the above notations the diffusion convection equation is dgi + d(ugi) = 1 d_ dt dx 3 dx

where

Pi+1

is a slowly varying function over much of momentum space, the simplest natural sub-grid model for fp) assumes q is constant inside Ap,; that is, we can apply a piecewise power-law model for fp) and use a logarithmic spacing in the momentum grid. Then, for example,

with obvious extension to gi and the other terms in the diffusion convection equation moment equations, where fi = f(Pi) = (+1/Pi)f+b di = Pi+1/Pi . (4.23.18)

Using the ratio gi/) one can derive the index qi once ni and gj are known.

Fig 4.23.1. shows results from an initial test of the scheme inside a TVD hydrodynamics code and, for comparison, the equivalent simulation with a conventional finite difference scheme for solving the diffusion convection equation (Gieseler et al., 2000). The simulations involve evolution of an initial Mach 10 shock discontinuity with an initially uniform CR pressure equal to twice the