[cyou]

I am now studying the immersed boundary methods
Do anyone know the ways to define the ghost cells?
Could you recommend me some papers?

Thx

[cyou]

nobody knows?...

[arjun]

Originally Posted by cyou

nobody knows?...

I definitely know as I have written immersed boundary based codes.

Before I answer your question, could you answer mine. Do you know what role ghost cells play in immersed boundary method. Do you understand the purpose of ghost cells?

[cyou]

Originally Posted by arjun

I definitely know as I have written immersed boundary based codes.

Before I answer your question, could you answer mine. Do you know what role ghost cells play in immersed boundary method. Do you understand the purpose of ghost cells?

To arjun

Thank you very much for your reply
Yes, I have read numbers of papers about the immersed boundary.
And I have did some simple simulation with the immersed boundary.
But the identification of the boundary and the interpolation are big problems for me now, especially the part of coding.

First, I want to define the marker points on the boundary
Then, find the ghost cells or the forcing cells
And then, find the image cells
Last, use the marker points, ghost cells(or forcing cells), and image cells to do the interpolation to add the body force

Is that right?

Thx

[arjun]

Originally Posted by cyou

To arjun

Thank you very much for your reply
Yes, I have read numbers of papers about the immersed boundary.
And I have did some simple simulation with the immersed boundary.
But the identification of the boundary and the interpolation are big problems for me now, especially the part of coding.

First, I want to define the maker points on the boundary
Then, find the ghost cells or the forcing cells
And then, find the image cells
Last, use the maker points, ghost cells(or forcing cells), and image cells to do the interpolation to add the body force

Is that right?

Thx

Yes and no.

Now at this moment, lets leave the concept of ghost cells aside. Assume that the solver does not use ghost cells, only uses immersed boundary method.

First of all understand that forcing is applied to momentum equation so that the regions of mesh where solid lies , the velocity is same as that of the solid.

This takes care of momentum equation. Second thing we have is continuity equation, that decides the pressure in the calculation domain.

This pressure equation is solved in whole domain , that is fluid and solid, just as it would be solved for normal solver. Nothing special here.

This shall work fine as a normal immersed boundary solver.


Now the question arises, why do we need ghost cell method at all , if that all could work fine?

The reason is that there are some errors remaining in momentum equation. Lets look at the errors that could be there.

When you apply descretization, you go through all the cells and calculate diffusional and convection terms based on u,v,w values of cell and neightbour cells.

what happens if neighbour cell happen to be a solid cell?? What happens with the solid-fluid interface lies somewhere between this fluid and solid cell??

In this situation the diffusional term for momentum equation would be wrong. Because diffusional term is equal to K(u_neigh - u_self) / d, where d is distance between cell centers.

If solid boundary lies between cell centers then correct diffusional term is K(u_neigh - u_self) / d_actual)

This d_actual < d , that means diffusional term as calculated by solver is lower than the exact value that shall be used.


Further it shall be noticed that (u_neigh - u_self)/d_actual is nothing but velocity gradient.

Now lets call

D_correct = K(u_neigh - u_self) / d_actual)
D_solver = K(u_neigh - u_self) / d

NOW HERE IS CRUX OF GHOST CELL METHOD:

Change u_neigh to u_neigh_ghost such that

D_solver = K(u_neigh_ghost - u_self) / d = D_correct

That is just change the value of velocity at cell center of solid cell to a value such that when solver calculates the gradient or diffusion term it comes out to be correct value.

Now lets come to your question: How to find ghost cells?


To find ghost cells loop through all the cells in mesh, if the cell is fluid cell, then check its neighours. If any of its neighbour is inside solid then this solid cell is ghost cell. Find the new velocity values for this ghost cell such that when gradients are calculated they are correct ones.

How to mark solids in whole domain is different question. It could be done by various methods. I won't mention how I am doing it because I am not allowed to discuss it.

[cyou]

Originally Posted by arjun

Yes and no.

Now at this moment, lets leave the concept of ghost cells aside. Assume that the solver does not use ghost cells, only uses immersed boundary method.

First of all understand that forcing is applied to momentum equation so that the regions of mesh where solid lies , the velocity is same as that of the solid.

This takes care of momentum equation. Second thing we have is continuity equation, that decides the pressure in the calculation domain.

This pressure equation is solved in whole domain , that is fluid and solid, just as it would be solved for normal solver. Nothing special here.

This shall work fine as a normal immersed boundary solver.


Now the question arises, why do we need ghost cell method at all , if that all could work fine?

The reason is that there are some errors remaining in momentum equation. Lets look at the errors that could be there.

When you apply descretization, you go through all the cells and calculate diffusional and convection terms based on u,v,w values of cell and neightbour cells.

what happens if neighbour cell happen to be a solid cell?? What happens with the solid-fluid interface lies somewhere between this fluid and solid cell??

In this situation the diffusional term for momentum equation would be wrong. Because diffusional term is equal to K(u_neigh - u_self) / d, where d is distance between cell centers.

If solid boundary lies between cell centers then correct diffusional term is K(u_neigh - u_self) / d_actual)

This d_actual < d , that means diffusional term as calculated by solver is lower than the exact value that shall be used.


Further it shall be noticed that (u_neigh - u_self)/d_actual is nothing but velocity gradient.

Now lets call

D_correct = K(u_neigh - u_self) / d_actual)
D_solver = K(u_neigh - u_self) / d

NOW HERE IS CRUX OF GHOST CELL METHOD:

Change u_neigh to u_neigh_ghost such that

D_solver = K(u_neigh_ghost - u_self) / d = D_correct

That is just change the value of velocity at cell center of solid cell to a value such that when solver calculates the gradient or diffusion term it comes out to be correct value.

Now lets come to your question: How to find ghost cells?


To find ghost cells loop through all the cells in mesh, if the cell is fluid cell, then check its neighours. If any of its neighbour is inside solid then this solid cell is ghost cell. Find the new velocity values for this ghost cell such that when gradients are calculated they are correct ones.

How to mark solids in whole domain is different question. It could be done by various methods. I won't mention how I am doing it because I am not allowed to discuss it.

Thank you very much for your detailed instruction.
I have come to understand something I have never thought before. Thank you very much...

It is not appropriate to ask you to introduce your methods in details now. But could you please recommend some useful papers or codes which can be shared on the internet?

[arjun]

Originally Posted by cyou

Thank you very much for your detailed instruction.
I have come to understand something I have never thought before. Thank you very much...

It is not appropriate to ask you to introduce your methods in details now. But could you please recommend some useful papers or codes which can be shared on the internet?

Search for Muzaferija S and Peric's paper for Ghost cells it covers in detail how it is implemented.

Or just search ghost cell method and read the paper while remembering the outline as I said. Now that you know the basic of it the papers will make sense to you.

[cyou]

Originally Posted by arjun

Search for Muzaferija S and Peric's paper for Ghost cells it covers in detail how it is implemented.

Or just search ghost cell method and read the paper while remembering the outline as I said. Now that you know the basic of it the papers will make sense to you.

Thx!!!
Hope everything goes well with your research!

[arjun]

i forgot to add, if you visit www.inavier.com you could download the inavier demo version that supports basic immersed boundary method on unstructured mesh.

I am going to improve it further in future when i get time, but for time being it is not bad.

[cyou]

Thx

And one more question: do you code the programm to do the simulation or make some modification on the commercial CFD programm? Which way do you recommend?

[cyou]

Originally Posted by arjun

i forgot to add, if you visit www.inavier.com you could download the inavier demo version that supports basic immersed boundary method on unstructured mesh.

I am going to improve it further in future when i get time, but for time being it is not bad.

Thx

And one more question: do you code the programm to do the simulation or make some modification on the commercial CFD programm? Which way do you recommend?

[arjun]

Originally Posted by cyou

Thx

And one more question: do you code the programm to do the simulation or make some modification on the commercial CFD programm? Which way do you recommend?

As far as I am concerned, I have written the in house code so I control it fully.

The code that I asked you to download is not the in house code, it is written as a part of hoby just for fun and learning. since inavier is also written by me, I can add models into it easily as I control the full source code.

Modifying commercial code is difficult but not impossible.

The major difficulty is the part where you mark the part of mesh as solid. I call it solid marking algorithm. It involves 3d range or nearest neighbor search. Because of this reason it could kill the efficiency of algorithm.

If you could do that by using external program and couple it with commercial code then modifying the commercial code is not difficult.


in iNavier the 3d search algorithm library is written by me, which is slightly different than other 3d search algorithm. normally people use kd - tree search for this purpose. I use Hyper-kd tree search algorithm that I developed. It is much faster than normal method.


My advise, do not hurry, sit down and think and think and think. weigh your options if you could really do it with commercial codes.

[cyou]

Originally Posted by arjun

As far as I am concerned, I have written the in house code so I control it fully.

The code that I asked you to download is not the in house code, it is written as a part of hoby just for fun and learning. since inavier is also written by me, I can add models into it easily as I control the full source code.

Modifying commercial code is difficult but not impossible.

The major difficulty is the part where you mark the part of mesh as solid. I call it solid marking algorithm. It involves 3d range or nearest neighbor search. Because of this reason it could kill the efficiency of algorithm.

If you could do that by using external program and couple it with commercial code then modifying the commercial code is not difficult.


in iNavier the 3d search algorithm library is written by me, which is slightly different than other 3d search algorithm. normally people use kd - tree search for this purpose. I use Hyper-kd tree search algorithm that I developed. It is much faster than normal method.


My advise, do not hurry, sit down and think and think and think. weigh your options if you could really do it with commercial codes.

Thank you very much for your suggestion!!!

[cyou]

Originally Posted by arjun

As far as I am concerned, I have written the in house code so I control it fully.

The code that I asked you to download is not the in house code, it is written as a part of hoby just for fun and learning. since inavier is also written by me, I can add models into it easily as I control the full source code.

Modifying commercial code is difficult but not impossible.

The major difficulty is the part where you mark the part of mesh as solid. I call it solid marking algorithm. It involves 3d range or nearest neighbor search. Because of this reason it could kill the efficiency of algorithm.

If you could do that by using external program and couple it with commercial code then modifying the commercial code is not difficult.


in iNavier the 3d search algorithm library is written by me, which is slightly different than other 3d search algorithm. normally people use kd - tree search for this purpose. I use Hyper-kd tree search algorithm that I developed. It is much faster than normal method.


My advise, do not hurry, sit down and think and think and think. weigh your options if you could really do it with commercial codes.

Long time no "see".
Hope you everything goes well

I have tried to write some codes with IB methods.
But there are few questions need your help.

1, If I simulate a 2D square whose boundary is coincide with mesh grid. Then there is no need to do any interpolation. Using IB method means define the velocity zero at the grids which are located in the square, after solving momentum equations and before solving the poisson equations. Am I right?

2, When the boundary is not coincide with the mesh grid. We need to do interpolation at the "ghost cell"( one layer grids between solid and liquid). But How about the grids in the model which means solid. Do we just define the velocity at the solid grid zero? Or we do not have to do anything to do with them?

Waiting for your kind reply

[calvinhobbes]

i am trying to write a solver based on ghost cell method proposed by Mittal 2008 for incompressible viscous flow usingg matlab for time being which will be converted in C later.i coupled the ghost nodes and fluid nodes while solving the matrix. In doing so, u- momentum and v- momentum pose no problem. but solving pressure poisson equation is becoming extremely dificult. No iterative solvers that comes with MATLAB can solve the matrix. I checked the conditional number of PPE which is around e^16. Any suggestion for this problem????And why this is happening?

[arjun]

Sorry long time and I did not see your post.

Originally Posted by cyou

Long time no "see".
Hope you everything goes well

I have tried to write some codes with IB methods.
But there are few questions need your help.

1, If I simulate a 2D square whose boundary is coincide with mesh grid. Then there is no need to do any interpolation. Using IB method means define the velocity zero at the grids which are located in the square, after solving momentum equations and before solving the poisson equations. Am I right?

Yes. Indeed you do not need any interpolation as you can close momentum eq. terms exactly.

Originally Posted by cyou

2, When the boundary is not coincide with the mesh grid. We need to do interpolation at the "ghost cell"( one layer grids between solid and liquid). But How about the grids in the model which means solid. Do we just define the velocity at the solid grid zero? Or we do not have to do anything to do with them?

Waiting for your kind reply

Understand this part that for any live or flow cell whose any neighbor is a solid cell, you can not close the momentum equation. If you are in a finite volume code then you will need convection term and diffusional term.

You can set the convectional term to zero (or really compute from interpolated value such that normal flux is zero). What really very important is that the diffusional term is set correctly.

This is what ghost cell method does. It sets the velocity in ghost cell such that when the diffusional term is calculated it gives you right value accounting for the surface of solid.

We noticed that this becomes very very important in flows with high Reynolds number with solids rotating. For example if this term is not correct your lift forces will be way off.

[arjun]

Originally Posted by calvinhobbes

i am trying to write a solver based on ghost cell method proposed by Mittal 2008 for incompressible viscous flow usingg matlab for time being which will be converted in C later.i coupled the ghost nodes and fluid nodes while solving the matrix. In doing so, u- momentum and v- momentum pose no problem. but solving pressure poisson equation is becoming extremely dificult. No iterative solvers that comes with MATLAB can solve the matrix. I checked the conditional number of PPE which is around e^16. Any suggestion for this problem????And why this is happening?

You matrix may be all Neurmann BC and thus singular. Check this part.

[calvinhobbes]

Originally Posted by arjun

You matrix may be all Neurmann BC and thus singular. Check this part.

But for fractional step method proposed by Von Neumann boundary condition is imposed on all boundaries. Does this mean the PPE may hav more than one solution? Also how to solve PPE in this kind of situation where Neumann boundary condition is applied?

[arjun]

Originally Posted by calvinhobbes

But for fractional step method proposed by Von Neumann boundary condition is imposed on all boundaries. Does this mean the PPE may hav more than one solution? Also how to solve PPE in this kind of situation where Neumann boundary condition is applied?

Yes indeed. You might have to check the documentation of matlab as to how this is done. I have used matlab exactly 2 times in life so do not know much about it.

[o_mars_2010]

Hello Mr. Arjun

I have developed a code to simulate a flow past a circular cylinder using Immersed Boundary Method.
I used the interpolation method by Balaras (2004).
However, the recirculation length L/D at Re=40 is around 2, however, from literature i found that it is in the range 2.2~2.4.
Could you please tell me what could be the reason??
(I used different grid size, however i got the same results)
In addition, i am using Tecplot to calculate the drag coefficient.
For example: If the Cylinder radius is for example R=0.05 & the grid size is dx=0.01
I could not get accurate results for the drag coefficient if the the calculation is taking place with the same radious.
But, i get accurate results for the drag coefficient if i used a radious smaller than the actual radious such as R=0.0485.
What could be the reason????