newbie question: finite difference, immersed boundary method
 

Dear forum members,

I would like to solve a PDE arising from a variational problem related to a mass conservation problem plus regularization terms. The preferred solution method is finite differences due to the simplicity of implementation. The application area essentially is the extraction of flow information from image sequences depicting fluid motion in curved tubular structures (vessels).

The question I would like to address is how to handle the domain boundaries within a FD formulation appropriately where the boundary is given by a 2D curve that is not necessarily aligned to the uniform grid.

So far I came across the Immersed Boundary method that as far as I understand introduces additional force terms to the PDE dependent on a set of arbitrarily placed points that describe the boundary.

My question is essentially if the IB method is the proper approach to this problem, or are there other possibilities that work with FD and my be an alternative? Additionally if anyone knew about an example C/MATLAB code for IB, that would certainly be of help.

In case, to make the problem simpler, and a boundary with edge aligned segments is considered is there a recommended way to handle the BCs? For rectangular boundaries the addition of "ghost points" (virtual domain points for which setting fixed values at each iteration depending on the values of neighboring cells at the previous step enable the handling of common BC types in an elegant manner) seems as a good solution, however for domains containing "step" boundary segments the approach cannot be directly applied, as a ghost cell may have multiple neighbours within the domain and thus its value cannot be set to satisfy multiple BCs at once. Is there an established way of using such ideas with non rectangular domains?

Thanks, and looking forward to your feedback,
Peter

Hi, the IB technique is only one possibility among others... but you can use finite volume method that is very simple and can be generalized properly over volumes of general shape without using complicate mapping between physical and computational domain. As suggestion, you can read the book of Peric and Ferziger.

If you have simple and smooth geometry a transfinite grid interpolation along with coordinate mapping works well with FD methods.

The book by Andersson describes this well and it is quite easy to implement.

If the geometry becomes more complicated then there is a high risk that you won't be able to use this method however. Here unstructured FV methods might be better, but they are, imho, more difficult to code.

gl hf

The purpose of Immersed Boundary methods as its name indicates it (immersed) is to handle obtsacle within the domain and not really the boundaries of the domain,even if it is always possible. But I think in the foundation paper of Peskin he designed it to handle moving artery wall for blood flow.
But anyway according to your problem as it has been suggested by Filippo and Fordprefect, herearethe different options:

1) If you really want to use FD you have to use coordinate transformation (mapping) to switch from physical domaintocomputational domain
2) use finite volume designed for non orthogonal.the book of Peric is a must and you will find some fortran code to deal with such configuration.
3) you can use also Finite Element method,but muchmorecomplicated to implement,so I would advise you to avoid this option.

As Filippo stated it above 2) is not necessarly more complicated than 1) but it is really more powerfull.
So for me 2) wins!!

Good luck

FD vs FV
 

Thanks a lot for the replies from all of you!

Actually the question has arisen in the context of implementing an image processing algorithm that is based on a variational PDE approach where the original paper does not address the question with nonrectangular domains, however for practical applications with vessels it has to be considered while the paper only focuses on cases with straight tube arrangements. The preference for FD comes from the fact that the paper uses a finite difference discretization.

The domain transformation approach mentioned in Chapter 5. of Andersson seems feasible to me, however to have a complete view of the possibilities I would like to ask if any of you knew a book / paper in which the handling of generic curved boundaries is addressed on the original FD domain by using more sophisticated handling of BCs. In 8.1.1 of the Peric book a method is briefly described where the boundary is approximated with a "stepwise" one that is aligned to the edges of the FD grid. Unfortunately no details are disclosed (only a reference is given to "Large-eddy simulation of turbulent boundary layer flow over a hemisphere" by M Manhart, H Wengle that I cannot access). Do you have any suggestion for such methods?

Thanks again for your attention.

you can find many papers (see JCP) that treat a curvilinear boundary with step-wise rectagles... but I do not suggest to use such a method, it is very low order accurate and requires an expensive clustering of cells...

As the IB method is concerned, I suggest these three papers:
1) M.D. de Tullio, A. Cristallo, E. Balaras, R. Verzicco, Direct numerical simulation of the pulsatile flow through an aortic bileaflet mechanical heart valve, Journal of Fluid Mechanics, Vol. 662, 2009, pp. 259-290 (ISSN 0022-1120).
2) M.D. de Tullio, P. De Palma, G. Iaccarino, G. Pascazio, M. Napolitano, An immersed boundary method for compressible flows using local grid refinement, Journal of Computational Physics, Vol. 225, 2007, pp. 2098-2117 (ISSN 0021-9991).
3) http://www.annualreviews.org/doi/abs....061903.175743



and if you want reading a method mixing IB and FV on unstructured grid:

http://onlinelibrary.wiley.com/doi/1...d.278/abstract

Good work ;)

Unfortunately I have my Ferziger/Peric book at home and can not check what you are referring to. Perhaps I misread you but do you want to implement stepwise boundaries? This is usually not preferred since your numerical surface normal is not aligned with your real surface normal. I guess it depends on what type of equations you are trying to solve however.

I do not understand well your question.
First you have to discretize your equations on the domain.
1) you choose to map your physical domain into a computational one where you will apply your FD method
2) you will manage your boundary conditions excatly how you use to do it.

or
1) you choose FV method and then you grid your physical domain with non orthogonal control volumes and you discretize your equations on these non orthogonal volumes (see Peric)
2) you manage your boundary conditions how it should be with finite volume (see Peric).

Avoid the use of stepwise boundaries. It is not how we should do CFD in 2013 ;)

The issue is very complicated.

The reason as to why there is not much available on it is that there is no perfect way of dealing with general boundaries in immersed boundary codes.

Having written immersed boundary solver for both cartesian grid and in general unstructured grids, I can tell you that you will not go far if you just follow what is mentioned in the books.

It took me almost 8 months of effort to get drag and lifts to match experimental data for our calculations. In this I ended up doing things the way NOT mentioned in any literature.

I won't be able to share what I exactly did to make it work because my last company is in process of patenting it. So until they are done with it, I can not tell you what it was.

But here are few pointers for you:

(1) If you have cartesian type mesh then try to find out extactly where in X,Y and Z axis the solid boundary lies.
(2) Write the various descretization operators while accounting for position of boundaries.
(3) if you can not do above two things then spend your time with Ghost cell methods and the idea behind it, that would give you the best results specially when things are moving and rotating.
(4) You can apply the forcing in Matrices too by just modifying the source terms. I did that compared to calculating forces and distributing them etc etc. Read (5) as to why.

(5) Think that the solid regions are the regions of the mesh where velocity field is known. Now how you force them is upto you. Ultimately what you are doing is that calculating the source terms for momentum equation such that when you solve them in matrix form the velocities that matrix solver would give is same as solid velocities in solid regions.

(6) The major headache and time consuming part is solid marking in immersed boundaries. In this if things are moving then you will have to mark the solid and fluid regions of the mesh after every movement. This could be very problematic in non uniform cartesian meshes and in unstructured meshes. This typically involves solving nearest neighbor problems which are time consuming. The situation is grave in explicit codes where time step is small and courant dependent.

(7) If used properly and immersed boundaries are very powerful tool.


On contrary, I intend to develop on stepwise boundaries in coming years. :-D , this might be the way to go.

I agree, no books or papers will give you the tricks to make an IB code properly working...

We can not advise someone to use IB just to deal with his boundaries which are curved !

The only justification of using IB compared to body fitted grids is if we have moving boundaries or moving obstacles within the domain. In other cases using FV on arbitrary shape control volumes is much more powerful and accurate.
Or at least use transform coordinate, mappe your domain and use FD if you want.
But using setpwise boundaries combined with IB is definitely not the best option to deal with a boundary. By the way don't you know that most of IB methods are just first order?
After this everyone can do what he wants in life ;)

Finite difference, Immersed boundary method
 

I am combining the Immersed Boundary Method and finite difference method. However, for the geometry corner, it is difficult to determine which velocity is used to interpolated, because u or v forcing are near the both two boundary line. Do you have the same experience about this problems?

I have met the similar problem. How do you solved that?

The way to handle such situations is basically based on using multiple values. In practice, the cell geometry is fixed but for a given variable you will need as many internal values as many fluid cells will need them.

Consider for example the corner cells for a square cylinder with different temperatures on different sides. You can also have a fully structured grid with the cylinder surface aligned with the grid but, in the end, each corner solid cell has to store information for each of the neighboring cells.

Ibm
 

But how do we use so much stored information? In other words, which measure do you take to handle the problem?

From the practical point of view the answer is: I don't know, as i never programmed it.

In general, if i would like to make it completely transparent to the solver i would add an integer index for each multi-valued ghost cell, each index pointing to the values for the given neighbor fluid cell. then i would treat such multiple values as different unknowns in the system.

Then you would just need an if clause (based on a multi-value flag and the integer index) to fill the system matrix.

Actually i don't know if the ghost cells are usually treated implicitly or not. If not then the previous procedure is even more trivial (i guess).

However, i'm not an expert on immersed boundary. An article treating such multi-value method is the following:

http://www.sciencedirect.com/science...45793007000394

which, by the way, might be different from my guess

Consider also the following paper:

http://www.sciencedirect.com/science...21999108000235

I think its description of the overall immersed boundary method is made crystal clear, and it also describes the case of thin surfaces. I might not appreciate their somehow arbitrary treatment of sharp-edges but, still, it seems to work.

Ibm
 

Thanks a lot. I will try to read these paper

Problems
 

For instance, rectangular with bottom line connected to the bottom boundary, the forcing in horizontal line is very small and that in vertical line is much larger.

Do you have similar experience ?

IBM problems
 

Hi,


For instance, rectangular with bottom line connected to the bottom boundary, the forcing in horizontal line is very small and that in vertical line is much larger.

Do you have similar experience ?

I admit that i have not exactly understood your question. However, i have no specific experience in IB programming, so i can't do anynithing more useful than suggesting some reference.

At present IB method is largely used and preferred over other methods.
IB method is very old, Peskin proposed it in 1972 (2d flow in the heart) and renewed interest has been gained in this method for the fact that simple cartesian CFD codes can be used for managing complex geometries.

See http://www4.ncsu.edu/~zhilin/TEACHIN...8Z/Peskin1.pdf

The method was modified in several ways (e.g. the immersed interface method by Leveque) and you can see a review at http://www.annualreviews.org/doi/abs....061903.175743

Other than using IB, curvilinear surfaces were sometimes introduced in the BC on cartesian structures grids by proper interpolation along grid-step approximating the surface