Robin type boundary conditions for a Fokker-Planck equation

Hi,

I wrote a solver for a Fokker-Planck equation a couple of month ago. The equation is a minimalistic model for unidirectional pedestrian flow in a corridor (rectangle) - people enter the corridor on the left and exit on the right. The upper and lower boundaries correspond to walls. We used an HDG discretization for the diffusive operator and a DG for the convection. The boundary conditions were a bit tricky back then and Joachim sent me a 'quick fix' (which worked fine )
Now I want to consider a more complicated geometry - my corridor has a bottleneck - and this quick fix is not working anymore (people enter at the vertical boundaries or even all boundaries,....) So in principle I need to figure out how to include the Robin boundary conditions correctly.

The equation is as follows (sorry for the latex, but I was not able to figure out how to type equations here). Given a potential (either b = (vmax,0) or the gradient of the function satisfying the Eikonal equation) and inflow and outflow rates fin and fout solve:

u_t = div(\nabla u + u (1-u) b)

with BC
Entrance (inflow): (\nabla u + u (1-u) b) n = fin (1-u)
Exit (outflow): (nabla u + u(1-u) b) n = fout u
Walls: (\nabla u + u (1-u) b) n = 0

The 'old fix' looks as follows:

I changed my geometry and labelled the boundary segments as 'inflow', 'outflow' and 'wall'. But I don't use this anywhere - which I guess is the reason why it is not working

I also attached the full code in case you want to have a closer look at it. In principle I define the geometry (corridor with bottleneck), solve the eikonal equation using a fast sweeping method proposed by Quian et al, SIAM J Num Anal 45(1) 83-107 2007 (I could do this more efficiently I assume but that's not the most important part now), take the gradient of the potential as velocity field and then solve the parabolic FPE using an IMEX scheme - fast_sweeping_eikonal.py is the main, while functions.py includes all functions used.

Thanks a lot
Marie-Therese

Hi,

you are right, you should use the boundary names.
This line restricts the BFI to the "outflow" and "inflow" boundaries. Thus you have the desired boundary condition on the "wall" as well.

Best,
Christoph

Dear Christoph,

thanks a lot - works perfectly

I have two more questions (sorry)
*) How can I include homogeneous Dirichlet bc at the outflow using the HDG-DG discretization. I know that this is not the best idea, but my colleague wants to see it. I would expect the formation of a boundary layer in this case. I tried

But then I get a rather weird behaviour at the exit. Is there anything else I have to set ? Maybe for the diffusion operator, or do I miss something ...

*) I would like to make a movie from the simulation (to convince my colleague in the US, that these boundary conditions are no good indeed ). How can I do that ?

Thanks a lot again - I really like the new NgSolve Interface. Great job !
Marie-Therese

Hi,

I'm not sure if such a boundary condition makes any sense at an outflow boundary.

For a general flux "b" I would recommend to split the boundary BFI. Otherwise the IfPos wouldn't do what you want.
To make a video it is the easiest to make snapshots from python.
In functions.py: and make the snapshot afer the redraw in "FPIMEX4": The "internal" has "visoptions" and "viewoptions" to set some GUI parameters. For example to turn of "autoscale" to get the same scaling for the video. Afterwards you can call "ffmpeg" to generate the video:
Best
Christoph