Solving simple electrostatic problem

I’m fairly new to PDE solving domain. Being overwhelmed with what julia offers GitHub - JuliaPDE/SurveyofPDEPackages: Survey of the packages of the Julia ecosystem for solving partial differential equations, I could decide on the best way to solve my simple 2d electrostatic problem with Dirichlet boundary conditions. I want to find voltage (or E field) inside a circular region with almost parallel plates. Voltage obeys Laplace’s equation with the boundary conditions drawn in the following picture.

A have looked at Ferrite.jl and it seemed like an overkill for my problem. I have looked at @ChrisRackauckas 's awesome tutorial Solving PDEs in Julia - Nextjournal, it helped but there were so many packages listed there.

My wishes are the following:

  1. Good accuracy, the plates would be barely tilted with respect to each other (maybe a fine mesh would be needed)
  2. No external program meshing, I’d like julia to do it (as the geometry is fairly simple)


I would recommend Gridap .
Gridap uses the finite element method. You basically has to write the equations in weak form almost exactly as you would by hand. Just checkout the examples and tutorials and you will find something useful.


You can use Gmsh.jl to mesh the geometry within Julia (with the help of the gmsh library under the hood) for use with Gridap.jl.

Is there a user-friendly documentation for this package?
This is not the first time we see Julia packages that refer to something else, that refer to …

1 Like

Are you sure you want to have the boundary conditions like this? In fact, you have three electrodes now with a given potential. Is that what you want? I suspect you may want to have just the two electrodes, each with a given potential, and the boundaries modelled as if they were far away from the electrodes.

To learn how to use the package, follow the tutorial The tutorial is basically a few examples of mesh generation. The basics are pretty straight forward, you just need to get the hang of it. But for a simple mesh, just specify points and edges and surfaces (defined by edge loops).

@Paulo_Jabardo, thanks and sorry, because as far as documentation goes, it looks utterly unsatisfying.

Maybe one needs to start here and drill down.

You are correct that the documentation is not very good or clear. And linking to won’t help much. One of the main issues with gmsh, IMHO, is that it looks like a graphical mesh generator but in reality it is a scripting language for describing geometry that has a GUI. It is kind of difficult to get started and the best way is to go through the tutorial scripts which are well documented in the comments.

I just checked the julia tutorial and remembered that the julia tutorial codes have almost no comments (sorry for that…) but it is line by line translation of the python code that has comments Side by side comparison is fairly simple.

The reference manual has all the tutorial codes in gmsh’s own language and they are well documented even though the translation to Julia has some quirks.

That is the way to get started when using gmsh.

1 Like

Yes I’m sure, real ground in a cylinder with voltage V=0. This is a fairly common setup for beampipes in accelerators

1 Like

Looking at is feels like an overkill for my needs. The problem is 2-dimensional, ideally i would just want to see the mesh via Plots.jl or something similar. I’m not even sure if FEM is the right tool to use. Could it be advantageous to just use finite difference methods with very basic grid -like mesh?

Julia native meshing would be perfect, as I want to change the geometry slightly and perhaps do some optimizations with regards how to make fields look like I want

how can you plot the mesh then? in JuliaPlots/Makie?

Of course you can solve this just as a simple PDE, but given the nature of the problem, the geometry and your accuracy requirements I wonder if it might not be better to use something like a boundary element method? This way hopefully you’d just have to discretize lines.


One of my favorite methods is to represent the electrical potential as a radial basis function expansion and to solve for the coefficients. This is mesh free, you can just plop the basis functions down on centers either in a grid or with a layout that matches the geometry in some way.

Remember to place some centers outside the domain, this will help in getting the boundary conditions right.

@dlakelan, any citations or links you can provide on the radial basis function expansion approach?

I’d second the recommendation for BEM. Has anyone extended BEAST to potential theory, or are there other good packages?

This sound very interesting, would you mind giving some more details. Maybe some introductory texts, Julia packages, anything else?

Can’t you just take the Helmholtz kernel with k=0?

This sound very interesting, would you mind giving some more details. Maybe some introductory texts, Julia packages, anything else?

Don’t ask me, I know next to nothing about it. I’d start with the wikipedia page then try to look up sources close to your scientific field. Probably, like all Green’s function methods, it has 30 different names depending on the community. Keywords that might be related (again I’m just relaying what I’ve heard mentioned in connection to the topic): single/double layer potential, methods of moments

1 Like

Even with BEM you still need to learn how to generate a (surface) mesh, although in 2d this is just lines. Really, there’s not much alternative to learning Gmsh or some similar package in the long run.

I agree that BEM or related surface integral equation methods (e.g. Nystrom) could be the most computationally efficient methods for this problem, but they are finicky to implement well and have a steep learning curve. I suspect it will be a lot quicker (in programmer time) just to bite the bullet and learn Gmsh and Gridap, unless BEAST already solves exactly the problem you want. The problem you are describing is not too computationally demanding.

1 Like

At one point I had a number of citations of papers I had actually read, but my zotero instance is acting up. so I’m going to just point you to a variety of results from google, from which you can branch out. In my experience the “ill conditioning” of the matrices is less of an issue when using compact basis functions or those that decay rapidly (like gaussians). But I admit to not having extensive examples. I think of the conditioning problems as occurring when you need a lot of cancellation between different basis functions, which occurs when the basis functions are near-constant across large spatial regions. It can also be useful to include in the basis a set of functions that represent the partial derivatives of the main basis functions.

The essence of solving PDEs in general is to create some family of functions that are general purpose approximators, and then to either evaluate some integral/weak form of the equations, or set the equations to be exact at some set of points (colocation methods). In the RBF expansion the RBF expansions are the family of functions, and typically the colocation method is used at a set of points. In the FEM method low order polynomials on patches of space knitted together are the family of functions, and either collocation at nodes, or some kind of gaussian quadrature weak form is used to get the equations to solve.

The advantage of the RBF approach is that you don’t need a mesh, and typically you get exponential convergence with respect to the number of basis functions, also it works without a regular geometry as compared to some spectral methods which work much better for say periodic boundary conditions on simple regions of space.

I’ve actually done this in Maxima to solve some heat equations on simple geometries like the one described here, with as few as 5-20 basis functions and gotten good results.