Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Solver for Poisson's equation in (radial,polar) coordinates. #245

Open
wants to merge 23 commits into
base: master
Choose a base branch
from

Conversation

mrhardman
Copy link
Collaborator

@mrhardman mrhardman commented Sep 6, 2024

This PR gives the option of solving

$$\frac{1}{r} \frac{\partial }{\partial r}\left( r \frac{\partial \phi}{\partial r}\right) + \frac{1}{r^2}\frac{\partial^2 \phi}{\partial \theta^2} = \rho(r,\theta)$$

for $\phi = \phi(r,\theta)$ using two different methods.

First, Fourier transforms for $\theta$ and Gauss-Legendre finite-element methods in 1D for the resulting system of independent ODEs. A Fourier spectral object (separate from Chebyshev) is added. Dirichlet (zero) boundary conditions are imposed on $\phi$ at the radial boundary. In the test, the radial domain runs from $(0,L]$.

Second, a 2D tensor-product Gauss-Legendre finite-element method, using the functions originally designed for the Fokker-Planck collision operator. Some base level functions have to be moved to effect this change. The sparse-matrix version of the assembly where the periodic boundary condition is imposed is currently broken, so for this PR we resort to the dense-matrix version of the operator creation (for now).

These features (and their extensions) may be useful if more complicated field solvers are eventually required in moment_kinetics, for example, for solving for $\phi$ from vorticity on closed field lines, or solving for the electromagnetic fields.

An automatic test is added which runs in serial (even if multiple cores are passed).

…ian in cylindrical polar coordinates is required for a particular dimension. Make vperp and r the two dimensions for which this is true. Note that this does not generalise well to other situations with more complicated Jacobian factors.
…or non-trivial functions of polar angle (with phi = 0 boundary conditions only).
…ndrical coordinates (radial,polar), as opposed to (radial,vertical) coordinates used for the velocity space solvers in the Fokker-Planck operator.
… to enable them to be reused for a different set of coordinates.
…ndre tensor product elements (not using Fourier transform).
…of periodic function -- make sure data for lower row assembly comes from correct location.
…matrix using ic_global and icp_global indices only (no icsc index). This assembly is slower for larger resolution due to memory usage, but can be done correctly with existing indexing functions. Update automatic tests to test the 2D Poisson solve as well as the 1D+Fourier decomposition version.
… These 1D ODE solvers are useful to show the basics of how the more complex 2D PDE solvers are constructed, and to show how different types of boundary conditions are imposed.
@mrhardman mrhardman added the enhancement New feature or request label Sep 19, 2024
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.

1 participant