[MC97]

I've to perform a multi-fluid internal flow simulation with a code which intrinsically assumes periodicity on a given direction (say z) on a cartesian grid. Nonetheless, the problem I'm trying to solve does not involve this kind of condition, in fact the fluids are running out from the computed domain through the two z extreme edges (one of them at sonic speed, the other at almost sonic speed) after being generated inside it from a given constant source. The simulation evolves the two fluids profiles up to a quasi-stationary state, when it is stopped.

Due to the code structure, I can't change the coordinate set (hence, I need to keep z as the direction on which the two fluids are getting out of the computed domain). I was then trying to figure out a way to solve this problem without removing the periodicity assumption, which would require a painful amount of coding. I thought of adding a penalized buffer layer at both z direction extremes on which I could impose the boundary conditions I need.

The problem is I must have Dirichlet BCs for velocities, and homogeneous Neumann BCs for fluids density at z extremes, and this would lead to a constant accumulation of the two fluids in this buffer layer, which will at least give some troubles to reach a quasi-stationary state for the simulation to end. I thought of adding a densities sink inside this buffer layer for this reason but I'm not sure about how this could impact the simulation inside the "real" computed domain (since one of the two fluids is not actually flowing at sonic speed).

Hence my questions are:

- how do you think the densities sink inside the buffer layer would influence the simulation inside the "real" computed domain?
- do you have any other idea to deal with this problem?

Thanks a lot in advance to anyone who will pay attention to this.

[sbaffini]

I think it really depends from how far away you can put those boundaries with the relative buffer layers.

Which, however, leaves on the table the main question, the suitability of the code you choose for the given task. The only scenario I can think of where recoding the bc would mean more effort is for a code that is spectral in the z direction. But then, multi-fluid, sonic speed... uhm, you're gonna hurt you so much with this, are you sure you wont to do it with this code?

My honest suggestion is, if your research/work is about this adaptation itself and/or you have time to work it out, then probably proceed and find a way to make it work. Otherwise, you need (in this order):

1) A lot of luck
2) All the literature ever produced on buffer layers
3) A lot of luck again

If 1 and 3 work, you either find exactly what you need in 2 or whatever you try yourself will just work at the first try.

If, instead, the code isn't spectral in z, honestly, just substitute the buffer layer with proper ghost cells implementing the boundary condition you want and just don't solve for them.

[MC97]

@sbaffini thank you for your feedback!

>> "My honest suggestion is, if your research/work is about this adaptation itself and/or you have time to work it out, then probably proceed and find a way to make it work."

Part of my project is indeed adapting this code to the situation I need to model. The code has other features (it is a plasma turbulence code) which motivate the interest in doing this.

>> "If, instead, the code isn't spectral in z, honestly, just substitute the buffer layer with proper ghost cells implementing the boundary condition you want and just don't solve for them."

This would actually imply breaking the periodicity assumption on z. Which is cool and of course physically sounded, but I'd like to figure out another way to do this because the code is parallelized and this would need a fairly amount of re-structuring of the code basis (which wasn't meant to be the aim of my project, as I'd originally need just to work at "higher level"). Of course I'll do this if no other option will be available.

[sbaffini]

I don't understand yet if the code is spectral in z.

But, IF IT ISN'T, this is not really different from a buffer layer... of course with few caveats here and there, but nothing really complicated if the code has a decent structure.

In this way you break the periodicity, yes, but if your flow actually exits from both ends in z, I suspect it is something you want to break indeed.

On the other end, if it is spectral in z, the buffer layer is probably your best/only shot. Maybe you can find a first batch of references on the topic of buffer layers in CHQZ, and move from there. It has been done, I'm just not sure about the multi-fluid/compressible scenario.

[MC97]

>> "I don't understand yet if the code is spectral in z."
Sorry. No it isn't.

Anyway, ok, I think I got what you mean. Thank you again for your advice .

[sbaffini]

Just as an example, assuming your scheme is finite volume with only nearest neighbor cells in the stencil (which is 2nd order), and a simple 1D scalar convection equation:



If your code is periodic in z, somewhere, it will do something like (Fortran notation):

phi(:,:,nz+1) = phi(:,:,1)
phi(:,:,0) = phi(:,:,nz)

You, then, will leave those untouched, and add, afterward, something like:

phi(:,:,nz+1) = phi(:,:,nz) !Neumann at nz
phi(:,:,0) = phi_in !Dirichlet at 1

If the code is implicit, this gets slightly more akward, because you need a second part that will do two things, remove the implicit contribution of the original periodic bc and add your own.

If your code has more cells in the stencil, I just expect the code to have more cells outside the domain to properly handle them. Now, the difficult part just becomes using the proper bc (as mine was just an example). But I expect your code to have the bc you want somewhere else, you just need to recode that into your own.

Now, of course, if we generalize this to whatever method your code uses to assign periodicity in z, this is just actually extending the code to a non periodic z bc that you will actually code. But, from my point of view:

1) In a well written code this should be no more complex than adding a buffer layer
2) This is actually simpler to code correctly (formally speaking, with respect to the underlying code) than coding correctly the buffer layer thing
3) Is guaranteed to work if the original bc you are coping works

The parallel thing might sound scary, but you will find that in the bc you are coping as well, and note that it is not any sort of magic. Any given cell needs, at some point, its neighbor ones to properly compute things with a stencil. The thing would have been more complex if you needed to add periodicity yourself. But as you are just implementing a local bc, you just need to be sure to do it in a place which is guaranteed to be not overwritten by some other routine.

Of course, if you instead have to investigate this extension AND the buffer layer way, there is nothing else to add. You need to do it. I would just consider that, when you do stuff in a compressible code, like adding a sink mass term that istantaneously cancels all the mass just entering a given zone, you are actually breaking the wave physics with something instantaneous. Sometimes this is just bad. So, when you have no reference, my suggestion is: you either know your numerics very well or just try to stick to the physics as much as you can. In this case, in lack of a better alternative, the physics would be the actual bc.

[MC97]

@sbaffini thank you for the really good answer

I will consider this, you've been really helpful.