[DR] Reduced Dimensionality Canonical Poisson Bracket Species

[johnson452]

Summary

Based on the branch: can_pb_reduced_dims. This was original branched from main, but should be merged into main and gk-g0-app branches.

This branch extends the Canonical Poisson Bracket solver to handle 1x2v, 1x3v, and 2x3v cases. (Ser p2, Hybrid p1, and tensor p1/p2). The requirement that a reduced dimensionality can be used is if the hamiltonian does not depend on the coordinate (which also implies the the metric does not depend on this coordinate too).

Since the hamiltonian depends on coordinate-momentum pairs, the code assumes q^i is paired with p_i. Therefore, the ignorable coordinate is presumed to start at the highest index. For instance:

Coordinate Momentum Pairs

Original can-pb (unchanged): (1x1v) : q = [q^0], p = [p_0, p_1] (2x2v) : q = [q^0, q^1], p = [p_0, p_1] (3x3v) : q = [q^0, q^1, q^2], p = [p_0, p_1, p_2]

Reduced dimensionality: (1x2v) : q = [q^0], p = [p_0, p_1] (1x3v) : q = [q^0], p = [p_0, p_1, p_2] (2x3v) : q = [q^0, q^1], p = [p_0, p_1, p_2]

These will be demonstrated through example regression tests below. Otherwise, this works exactly the same as the Vlasov species, and along as you hand the correct cdim and vdim, the kernels for ignorable coordinates with can-pb will be automatically selected.

Tests

Example regression tests [Already working]

• 1D Sodshock (x, y) (1x2v) - y-ignored Flat Space regression/rt_can_pb_neut_bgk_sodshock_im_1x2v_p1.c
• Annular disk (r, theta) (1x2v) - Theta-ignored regression/rt_can_pb_bgk_surf_annulus_sodshock_im_1x2v.c

Planned tests [WIP]:

• Sphere (theta, phi) without perturbations only doing theta resolution and compare to full case (1x2v) - phi - ignored
• Toroidal (r, theta, phi) (2x3v) - phi - ignored
• Toroidal (r, theta, phi) (1x3v) - theta, phi - ignored
• Cylindrical (r, theta, z) (1x3v) - theta, z - ignored

Notes

• 2x3v p2 tensor kernels are too large, so they are not available.
• Tensor p2 will not be testable until merged into gk-g0-app.
• Important for connecting the can-pb neutrals to other species is respecting the coordinate-momentum pairs. For instance the configuration q = [q^0, q^2], p = [p_0, p_1, p_2] is not valid.

[JunoRavin]

One thing we should strongly consider at this point is defining a set of kernels for flat spacetime and the set of kernels for general coordinates. I imagine the flat spacetime kernels will be much, much cheaper because then the Hamiltonian is only a function of velocity space (and I think this is the right nomenclature to use; a set of kernels where H is only a function of v, and a set of kernels where H is a function of x and v).

[johnson452]

Without any spatial curvature there will be no constraint forces. Also, splitting E and B off also means no need for updating fluxes across the velocity space surfaces. To update existing Vlasov to rely on the Hamiltonian my thoughts for followup steps are:

Taking the existing bracket kernels and doing this requires modifying:

Maxima

• Deleting the updates for surfaces in vx, vy, vz (this is just a few lines of code)
• Modifying the bracket to assume \pfrac{H}{x^i} is zero (this is just a few lines of code too)
• Change H to expand over velocity space basis set only.

G0

• Copy the can-pb infrastructure to overwrite the Vlasov-advection. Remove h_ij_inv and det_h, resize surf_alpha, remove velocity space surface kernels calls, and project H onto a single cell in conf space over the velocity space.

With these changes, almost the entirety of the existing can-pb infrastructure should just work right away.

[JunoRavin]

Depending on how we want to approach the G0 side, one could generate unit h_ij_inv and unit det_h, so that then infrastructure such as project Maxwellian doesn't have to change.

Likewise, surf_alpha can stay the same size because the configuration space updates are all that are used and those are the first indices always (and then long term the surf_alpha array's velocity space components can be filled with whatever forces are present).

[ammarhakim]

I want to understand @JunoRavin's suggestions here better. I think the distinction need not be made and we can in general assume H = H(x,v). There are very few cases in which H = H(v) only: flat-space no mappings. We have kernels for those already.

[JunoRavin]

@ammarhakim My suggestion is related to future work for the Vlasov infrastructure expecting the Hamiltonian to be specified, always, and then we can utilize a single set of kernels for things like non-relativistic Vlasov vs. relativistic Vlasov (since the only thing that changes is H = 1/2 m v^2 vs. H = m c^2 gamma) and we can also look towards a future in which all the Vlasov infrastructure works with general geometry

[johnson452]

Update on regression test: Most new regression tests qualitatively produce the right answer for 1x2v, 1x3v, and 2x3v for hyb/tensor/ser. (see attached slides) However, very small errors (currently the order of 1e-12) of total density and total number of particles are creeping by the end of the simulation. Tomorrow I plan to narrow this down to either the collisional accuracy or issues with the Hamiltonian continuity in these new cases. Can_pb_RT_Test_Results.pdf