Skip to content

Implicit, Fokker-Planck, self collisions for ion species (squashed) - #383

Merged
johnomotani merged 3 commits into
masterfrom
implicit-FP-collisions
Jul 13, 2025
Merged

Implicit, Fokker-Planck, self collisions for ion species (squashed)#383
johnomotani merged 3 commits into
masterfrom
implicit-FP-collisions

Conversation

@johnomotani

Copy link
Copy Markdown
Collaborator

Replaces #353 - squashed to clean up git history.

The aim of this PR is to provide the ability to advance the equation

$$\frac{\partial F_i}{\partial t} = C[F_i,F_i] \quad (1)$$

with the implicit, backward-Euler timestepping method, i.e., solve

$$F^{n+1}_i - F^n_i = \Delta t C[F^{n+1}_i, F^{n+1}_i] \quad (2)$$

using the newton_solve!() functionality provided by @johnomotani. This feature is integrated into the main time advance functions of moment_kinetics so that the implicit implementation of the collision operator can be used with one of the IMEX time advance schemes also provided by @johnomotani.

An example input file that uses an IMEX advance to solve the equation (1) above is provided in examples/fokker-planck/fokker-planck-relaxation-implicit.toml. Note the addition of the input section [fokker_planck_collisions_nonlinear_solver] which provides the flags for newton_solve!(), and note the presence of a new @enum for selecting the type of ion advance, e.g.,

    [fokker_planck_collisions]
    use_fokker_planck = true
    nuii = 1.0
    frequency_option = "manual"

    [fokker_planck_collisions_nonlinear_solver]
    atol = 1.0e-10
    rtol = 1.0e-10
    
    [timestepping]
    nstep = 10
    dt = 1.0e-3
    nwrite = 1
    nwrite_dfns = 1
    kinetic_ion_solver = "implicit_ion_fp_collisions"
    type = "PareschiRusso3(4,3,3)"

A separate set of inputs is provided for the ion collision non-linear solve as for the electron non-linear solve in case different default inputs are necessary. At this time I have not tested implicit ion collisions alongside implicit kinetic electrons, but it is plausible that different defaults could be convenient for fast simulations. The list of possible kinetic_ion_solver options are given here https://github.com/mabarnes/moment_kinetics/blob/11e4f651bce559643aa2f179efb630fdf1c3cec2/moment_kinetics/src/input_structs.jl#L51-L65.

A series of test scripts have been developed in the building of this feature. In test_scripts/ImplicitCollisionsTest.jl we have the test scripts test_implicit_collisions() https://github.com/mabarnes/moment_kinetics/blob/11e4f651bce559643aa2f179efb630fdf1c3cec2/test_scripts/ImplicitCollisionsTest.jl#L63-L75 and test_implicit_collisions_wrapper() https://github.com/mabarnes/moment_kinetics/blob/11e4f651bce559643aa2f179efb630fdf1c3cec2/test_scripts/ImplicitCollisionsTest.jl#L217-L230, which can be used to solve equation one with backward-Euler and the initial condition

$$F_i(t=0) \propto exp\left(-\frac{(v_\| - v_{\| 0})^2 + (v_\perp - v_{\perp 0})^2}{v_{th0}^2}\right) \quad (3)$$

An example function call to step with one timestep of backward-Euler would be

julia> using Revise; using moment_kinetics
julia> include("test_scripts/ImplicitCollisionsTest.jl")
julia> test_implicit_collisions_wrapper(test_particle_preconditioner=true,test_numerical_conserving_terms=true,
    vth0=0.5,vperp0=1.0,vpa0=1.0, nelement_vpa=32,nelement_vperp=16,Lvpa=8.0,Lvperp=4.0, bc_vpa="none", bc_vperp="none",
     ntime=1, delta_t = 1.0, ngrid=5)

I have also extended the coverage of the automatic test scripts moment_kinetics/test/fokker_planck_tests.jl and moment_kinetics/test/fokker_planck_time_evolution_tests.jl. The changes in moment_kinetics/test/fokker_planck_tests.jl test

  • the preconditioner based on the backward-Euler solution of the equation
$$\frac{\partial F_i}{\partial t} = C[F_i,F_0] \quad (1)$$

with $F_0$ a specified distribution function, usually chosen to be $F_i$ at the last timestep;

  • and test the lowest level nonlinear backward Euler solve of equation (1) via equation (2).

The changes in moment_kinetics/test/fokker_planck_time_evolution_tests.jl

julia> using Revise; using moment_kinetics
julia> include("moment_kinetics/test/fokker_planck_time_evolution_tests.jl")
julia> FokkerPlanckTimeEvolutionTests.print_output_data_for_test_update("runs/path_to_your_run/path_to_your_run", write_grid=true, write_pdf=true)

In summary, this PR

  • adds the option for implicit backward-Euler advance with nonlinear ion self collisions,
  • provides a preconditioner that is developed from the linearised test-particle collision operator that might be used to make cheaper implicit model operators,
  • gives the option to use near-exact numerical error corrections to preserve the moments of $F_i$ (the default behaviour),
  • extends newton_solve!() to support iterations in the anyv region layout,
  • tests the feature with automatic tests,
  • and updates the tests to use the initial condition (3) to give a clear comparison in tests to published examples

Before merging, I still need to add

  • An error exception to catch use of vperp.bc="zero-impose-regularity", which is not supported.
  • Revert changes to set_defaults_and_check_values!() and make [fokker_planck_collisions.nonlinear_solver] a standalone namelist which is read directly from the input TOML in setup_fp_nl_solve(). The function setup_fp_nl_solve() should contain any Fokker-Planck specific defaults for the nonlinear solver. Remove the nonlinear_solver member from fkpl_collisions_input.
  • Include nl_solver_params.ion_fp_collisions in suitable timestep checks.
  • Set defaults to use preconditioner (currently preconditioner is not used apart from in test scripts) and supply sensible default way of choosing to update it (note that this requires some method of storing the preconditioner at all spatial points).

Before merging this PR I would request feedback on

  • structure of inputs to select the kinetic ion time advance options and the non-linear solver options
  • suitability of the tests
  • whether or not sufficient documentation is provided
  • whether checks with simulations with vpa, vperp, z advection are required
  • anything else that comes up in discussion.

mrhardman and others added 3 commits July 12, 2025 21:45
It should do this as the 'anyv' subblocks are subblocks of the
shared-memory block that is synchronized by _block_synchronize(), so
every subblock is synchronized. Done by calling
_anyv_subblock_synchronize() - although this is mildly inefficient
(MPI.Barrier() gets called twice), it is only used for debugging so this
should be OK.
Use `comm_anyv_subblock[]` instead of the incorrect `comm_block[]`.
@johnomotani johnomotani added the enhancement New feature or request label Jul 13, 2025
@johnomotani
johnomotani merged commit 04b42e1 into master Jul 13, 2025
21 checks passed
@johnomotani
johnomotani deleted the implicit-FP-collisions branch July 13, 2025 15:51
johnomotani added a commit that referenced this pull request Jul 13, 2025
This came from #327 originally, which turned out not to be needed.
johnomotani added a commit that referenced this pull request Jul 14, 2025
Revert removal of NCDatasets from CI, accidentally merged in #383
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancement New feature or request

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants