[time_evolving_new.zip](https://github.com/mabarnes/moment_kinetics/files/13459088/time_evolving_new.zip)

[mrhardman]

          [time_evolving_new.zip](https://github.com/mabarnes/moment_kinetics/files/13459088/time_evolving_new.zip)

time_evolving_1D2V_new.zip
Uploads including example input files and output .gif/.pdfs for a 1D2V sheath problem with C[F,F] and wall boundaries, and the corresponding Maxwellian relaxation problem (0D2V). These simulations suggest that commit c3b6261 behaves acceptably for 1) small enough time steps and 2) simulations where the stagnation point z=0 are eliminated with the choice of grids.

Originally posted by @mrhardman in #149 (comment)

[mrhardman]

Development in #149 since c3b6261 seems to have broken the long-time stability of sheath simulations with either numerical dissipation or fokker-planck collision operator. Examples to follow.

[mrhardman]

diss_test.zip
For example, switching to numerical dissipation for a speedy simulation, at long times we have numerical spikes and a crash of the simulation due to negative temperatures. The same qualitative behaviour is seen for nuii = 1.0

[mrhardman]

time_evolving_0D2V_diss.zip
An obvious bug is found if we run with only numerical dissipation in vperp and vpa. No sources are specified, but it looks like the solution is growing! @johnomotani Have I made a silly error in the input?

[mrhardman]

It appears that the "chebyshev_pseudospectral" option for vpa and vperp gives plausible behaviour but the "gausslegendre_pseudospectral" option gives anti-diffusive growth. @johnomotani What have we changed there recently?

This .gif is made by using the input above and changing the discretization options for vpa and vperp.

[johnomotani]

I suspect the problem is that the 'boundary terms' here https://github.com/mabarnes/moment_kinetics/pull/149/files#diff-6af4b87e0934f1c86e512a7fbaf16a705fcf73131def91f308817abbf9f14f24R1016-R1017 have the wrong sign - at least if I flip the signs of those then the evolution looks pretty much like the Chebyshev-discretization one. I haven't checked what those signs should be from the weak form discretization. If those signs are wrong, then I guess other boundary terms likely are as well?

[mrhardman]

I suspect the problem is that the 'boundary terms' here https://github.com/mabarnes/moment_kinetics/pull/149/files#diff-6af4b87e0934f1c86e512a7fbaf16a705fcf73131def91f308817abbf9f14f24R1016-R1017 have the wrong sign - at least if I flip the signs of those then the evolution looks pretty much like the Chebyshev-discretization one. I haven't checked what those signs should be from the weak form discretization. If those signs are wrong, then I guess other boundary terms likely are as well?

That sounds like a problem. Better tests for the second derivatives would be needed. Can you give a line number? The link doesn't resolve to specific file. If this is the issue, I don't think it affects the collision operator, which would explain why the Fokker-Planck relaxation test still passes. That means there is at least one other important bug.

EDIT: Thinking clearly, I did make tests of the boundary terms in https://github.com/mabarnes/moment_kinetics/blob/merge_fkpl_collisions/test_scripts/GaussLobattoLegendre_test.jl. This should mean that the extreme boundary terms are correct.

[mrhardman]

Improvements to the explicit boundary terms for the gausslegendre_pseudospectral option are here: d306695. This seems to solve the problem with numerical dissipation in isolation, but it does not seem to explain why sheath simulations are broken.

[mrhardman]

Using 6dfd0d9 we still see problems in distributed memory simulations, even with a simple diffusion operator in a 1D1V simulation. The .gif here shows how the distribution function calculated incorrectly. The .zip contains an example of a local simulation that appears to perform as would be expected.

diss_test.zip

The currently parallelisation of C[F,F] requires us to use distributed memory for z for a simulation to run in a time that is affordable. The fact that the collision operator tests pass suggests that the remaining bug is in the distributed memory part of the code, which is not tested at check-in.