CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > General Forums > Main CFD Forum

Matlab code for pipe flow

Register Blogs Community New Posts Updated Threads Search

Like Tree6Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 8, 2011, 05:40
Default
  #41
New Member
 
Vincent
Join Date: Jul 2011
Posts: 29
Rep Power: 15
VincentD is on a distinguished road
I was wondering prasanthnitt, did this solution work out for you?

Next it might be nice to change this code to 3D using multigrid method. Anyone up for that?

Regards,

Vincent
VincentD is offline   Reply With Quote

Old   October 8, 2011, 06:27
Default
  #42
Member
 
Prasanth P
Join Date: May 2009
Posts: 30
Rep Power: 17
prasanthnitt is on a distinguished road
Hi Vincent,

I am extremely sorry for not replying, was stuck with something and lost track of the post. I m still looking for better outflow boundary conditions. P=0 and dp/dn=-delta/dt both are working for me. Working in the sense the Poisson is converging. I have noticed one issue with the direct Poisson solver, for inflow outflow problems if we set dp/dn=0 at all the boundaries w still get some results out of it but an iterative procedure does not give a result for the same. I dono the reason for that.
I m yet to implement Multigrid methods for the 2d case. But it would be a good idea to extend the code to 3d. I will be doing it in a while.
I am still workin on the boundary condition issue. I need to find some for a turbulent case. People have used Outflow boundary conditions like convective boundary condition. Do any of you have any idea regarding that.

Cheers
Prasanth
prasanthnitt is offline   Reply With Quote

Old   October 18, 2011, 13:14
Default I need Matlab Code
  #43
New Member
 
Mohamad Hosein
Join Date: Aug 2011
Posts: 9
Rep Power: 15
Mh.R is on a distinguished road
hi

I need matlab code that solve NS on triangular unstructured grid finite volume method.

Thanks
Mh.R is offline   Reply With Quote

Old   October 18, 2011, 15:07
Default
  #44
Member
 
Prasanth P
Join Date: May 2009
Posts: 30
Rep Power: 17
prasanthnitt is on a distinguished road
Hi Mh. R,

I am also looking for one such code. In any language its fine not necessarily in matlab. Do u have an unstructured FV code in Fortran or C or C++? If yes please share it.

Cheers
Prasanth
prasanthnitt is offline   Reply With Quote

Old   October 19, 2011, 04:41
Default
  #45
New Member
 
Mohamad Hosein
Join Date: Aug 2011
Posts: 9
Rep Power: 15
Mh.R is on a distinguished road
Hi prasanthnitt

I have just fortran code that solve NS in structured grid, that written by my self, if you want it, I can send it for you.
Mh.R is offline   Reply With Quote

Old   October 19, 2011, 06:21
Default
  #46
New Member
 
Vincent
Join Date: Jul 2011
Posts: 29
Rep Power: 15
VincentD is on a distinguished road
Why do you want to use an unstructured grid?

Another option might be to use a structured grid with an immersed boundary method.

Regards,

Vincent
VincentD is offline   Reply With Quote

Old   October 19, 2011, 08:46
Default
  #47
New Member
 
Mohamad Hosein
Join Date: Aug 2011
Posts: 9
Rep Power: 15
Mh.R is on a distinguished road
Hi Vincent

My thesis is defined in this way.

Thanks for you suggestion and attention.
Mh.R is offline   Reply With Quote

Old   October 19, 2011, 09:46
Default
  #48
Member
 
Prasanth P
Join Date: May 2009
Posts: 30
Rep Power: 17
prasanthnitt is on a distinguished road
Hi Mh. R,

I have a structured code too, written by myself. My formulation is for an unstructured grid but for some strange reasons the unstructured results tat I get are not correct except for 2d-cavity flows. So I am working on that issue. Do let me know if you locate an unstructured incomp NS solver.
Wat problem are you trying to handle with the unstructured code?
Can you give the details?

Cheers
Prasanth
prasanthnitt is offline   Reply With Quote

Old   October 19, 2011, 15:11
Default
  #49
New Member
 
Mohamad Hosein
Join Date: Aug 2011
Posts: 9
Rep Power: 15
Mh.R is on a distinguished road
hi Prasanth

I want solving NS with sediment transport.
Mh.R is offline   Reply With Quote

Old   October 25, 2011, 07:26
Unhappy i am new
  #50
New Member
 
David LOBON
Join Date: Oct 2011
Posts: 1
Rep Power: 0
fisicas is on a distinguished road
i have read the post, i am new in the topic but my question is that i want to do a dynamic state of flow inside pipe this pipe a flux external (heat transfer) and need to do 2d in the pipe and time, do you you help me ? because the scrip Benjamin is for Solves the incompressible Navier-Stokes equations in a rectangular domain...but i need inside the pipe.... i apologize but i hope to help in the foro in two week

thanks
fisicas is offline   Reply With Quote

Old   November 12, 2011, 07:12
Default
  #51
New Member
 
Vincent
Join Date: Jul 2011
Posts: 29
Rep Power: 15
VincentD is on a distinguished road
Hi fisicas,

there is no such thing as a 2D pipe, a pipe is a 3D structure. I guess 2 parallel plates is about as good as it gets as for making a 2D pipe.

If your going to a 3D simulation of a pipe I recommend starting from the navier stokes equations in cylindrical coordinates. You can still treat the terms in a similar fashion as was done by Seibold.

I recommend leaving the heat equation out of it until you have build and verified your 3D pipe flow code. Adding the heat equation should not be a problem however.

Good luck!
VincentD is offline   Reply With Quote

Old   November 13, 2011, 10:46
Default
  #52
Member
 
Prasanth P
Join Date: May 2009
Posts: 30
Rep Power: 17
prasanthnitt is on a distinguished road
Hi Vincent,

Laminar Pipe flow is two dimensional.

Cheers
PP
prasanthnitt is offline   Reply With Quote

Old   December 22, 2011, 13:17
Default
  #53
Member
 
Peter
Join Date: Oct 2011
Posts: 52
Rep Power: 15
new_at_this is on a distinguished road
Hi,

I have modified the same MIT18086 code but I am having some issues so I was hoping that someone in this thread who is familiar with the code might be able to help.

The primary modification I did was to change the boundary condition at the lower boundary so that the velocity was equal to the upper boundary. This produces a symmetric flow with the line of symmetry at y = 0.5. Now I am trying to calculate the same flow, but only over half the domain by imposing the symmetry condition at the bottom boundary. It seems pretty strait forward but for some reason I cannot get the same answer. Does anyone here know how to do this?
new_at_this is offline   Reply With Quote

Old   May 26, 2012, 09:29
Default
  #54
Member
 
Michael Moor
Join Date: May 2012
Location: Ireland
Posts: 30
Rep Power: 14
michaelmoor.aero is on a distinguished road
Hi Everyone, I am about 80% (i hope) through writing m own 2d cfd code in matlab, using a backward staggered grid, hybrid differencing, steady incompressible, and Gauss-Seidel for solving the system of equations and for convection-diffusion flow between two stationary plates... I am aving trouble however in my first iteration, with the pressure correction system of equations not converging ( it is not strictly diagonally dominant)... my question is this:

Is it ok to allow the pressure corection not to converge in the first few iterations? will it get progressively more diagonally dominant as the solution develops? and hence hopefully allowing my overall solution to converge? or is this due to a poor choice of boundary conditions?

NI = number of nodes in x, NJ = number of node in y
i = velocity line, j = velocity line
I = scalar line, J = scalar line

I have attempted to use inflow outlfow bcs, the code asks the user for a stagnation pressure, which i then turn into a static one using bernoulli, and assign this to node I=1, (first scalar/pressure node), then I ask for a velocity input uin and vin, which i assign to i =2 and I=1 respectively, noting that u-cells solve from i = 3 onwards to NI-1 in the x-direction, and J=2 to J=NJ-1 in the y-direction, and v- cells solve from I=2 to I=NI-1in x, and j=3 to j=NJ-1 in y...

the walls have a no-slip bc, i.e u=v=0 along all nodes in the x direction, and nodes j=2 and j=NJ-1

the outlet then, located at i=NI. says that vNI-1 = vNI, uNI=uNI-1 *(sum uin/sum uNI-1) i.e multiplying by mass out over mass in, area and density are the same. sum inflow/sum outflow.

now versteeg says that the pressure gradient does NOT equal zero in the flow direction... so I have set the pressure at I= NI to zero...

the pressure correction is then solved from I=2 to I= NI-1 in x, and from J=2 to J=NJ-1 in y... and all my coefficient matrices have no zero terms along the diagonal... although they are not strictly diagonally dominant...

so again, if you would like to see the code i would love to put it up somewhere ( where I'm not sure as it has function files which might make it complicated to read if i put it in piece by piece here), and my questions are finally:

1. Are my boundary conditions correct? I have cut links in the momentum equations to boundaries and included the flux term from the inlet into the source term on the rhs of [A][x]=[b], and allowed uw=uw* in the pressure correction equation at the inlet, and then made no other adjustments to the pressure correction equations

2. If what I have done so far is correct, is it acceptable to use the unconverged values of the pressure correction in the solution and allow it to run it's course, or does the TDMA ( Tri - Diagonal Matrix Algorithm) solve for not strictly diagonally dominant matrices?

Apologies for all the info, I hope that it is clear and you smart chaps can help me out!!

Best Regards,
Michael
2.
michaelmoor.aero is offline   Reply With Quote

Old   October 2, 2012, 10:16
Exclamation Special treatment of pressure correction based on continuity conservation in a pressu
  #55
Member
 
Join Date: Mar 2009
Location: Istanbul, Turkiye
Posts: 47
Rep Power: 17
gemini is on a distinguished road
Hi everyone

Please look this paper for the implementation of pressure correction equation boundary conditions. It is very important for algorithms you are trying to develop.

SPECIAL TREATMENT OF PRESSURE CORRECTION BASED ON CONTINUITY CONSERVATION IN A PRESSURE-BASED ALGORITHM
http://dx.doi.org/10.1080/10407790190053842

Nice work to you
gemini is offline   Reply With Quote

Old   October 3, 2012, 07:11
Default
  #56
Senior Member
 
Join Date: Aug 2011
Posts: 272
Rep Power: 16
leflix is on a distinguished road
Hi guys,
I just followed your discussion and here are my comments.
I guess as you pointed out that the main reason your codes do not work is due to the way you implement your BC and also what kind of BC you use especially for the pressure.

Quote:
Originally Posted by prasanthnitt View Post

I have tried both P=0 and dp/dn =0. For the former case the poisson solver takes several thousand some times even lakhs depending on the domain size. For the latter case the poisson doesnt converge and doesn't diverge either.
dp/dn =0 is not a correct boundary condition when the velocity boundary condition is not dirichlet. But, I dono what is the issue with my solver.

Cheers
PP

dp/dn=0 is valid even if you do not use Dirichlet BC. Just consider that you use the same BC for U* (predictor step) than for U(n+1) corrector step. Then it leads also to dp/dn=0

For incompressible flows based on fractional step method (projection method, Chorin method) or SIMPLE like methods. dp/dn=0 is the right BC you should use for your poisson solver. Sometime extrapolation of pressure on the BC gives the same thing.
I know that some people claim that P=0 should be the right BC but for incompressible flows with such pressure-velocity coupling algorithm I'm a bit doubtfull. Anyway dp/dn=0 does work !
However as everyone know and as Vincent pointed out, a Poisson solver with such BC on all boundaries leads to a singular matrix.
Then you should use some specific recipes which consist in specifying one pressure node...a lot of posts on this site are concerned by this aspect..It is also not straightforward how to procede...
However some linear system solvers work well even with dp/dn=0 on all boundaries (sip solver, Jacobi, Gauss-Seidel, SOR)
For the later you just have to loop this way:

do while (res> epsilon)

loop on J
P(0,J)=P(1,J) here is included dp/dn=0
P(NI,J)=P(NI-1,J)
loop on I
P(I,0)=P(I,1)
P(I,NJ)=P(I,NJ-1)

do J=1,NJ-1
do I=1,NI-1

P(i,j)=......

end do
end do

end do

Some solvers do not admit such BC (dp/dn=0) and then they fail at least they do not converge.
Then change your solver. You can either do some stuffs in the Krylov subspace but it's a bit complicate and beyond my skills.


If you do not use FV method (but sometime it is also needed) you should correct your outlet velocity in order that the mass flow rate between inlet and outled is identical. It's a very common procedure which has been already mentioned in this discussion...


So to conclude, use dp/dn=0 with SOR or Gauss Seidel as indicated below, with velocity correction based on mass conservation and it should work.
leflix is offline   Reply With Quote

Old   April 3, 2015, 07:57
Default FVM using SIMPLE (flow through parallel plates)
  #57
New Member
 
deewakar
Join Date: Dec 2011
Posts: 18
Rep Power: 15
alligngr8 is on a distinguished road
Hello all,

Is this thread still active. I need help regarding finite volume method based code on Staggered Grid (SIMPLE method) for flow through parallel plates.

I need some guidance regarding boundary conditions.

Looking for reply

thanks
alligngr8 is offline   Reply With Quote

Old   April 3, 2015, 08:01
Default
  #58
Member
 
Prasanth P
Join Date: May 2009
Posts: 30
Rep Power: 17
prasanthnitt is on a distinguished road
Post your query
prasanthnitt is offline   Reply With Quote

Old   April 3, 2015, 09:53
Default here is my query
  #59
New Member
 
deewakar
Join Date: Dec 2011
Posts: 18
Rep Power: 15
alligngr8 is on a distinguished road
Hi
Thanks for reply.

I am facing the problem with boundary conditions. I will briefly describe my problem formulation followed by query.

I have a 2D domain. nx and ny are number of nodes along x and y direction. The approach is FVM and grid is staggered backwards.

U (x velocity nodes) is defined at nodes 1,nx along x and 1,ny+1 along y.
V (y velocity nodes) is defined at nodes 1,nx+1 along x and 1,ny along y.

The pressure nodes are nx+1 and ny+1 along x and y respectively.

(Please find attached snap of the same).
The equations being solved are
Quote:
a(J,i)*u(J,i)= a(nb)u(nb) -(P(J,I)-P(J,I-1)) .
% here nb means neighboring nodes.
% here I have written a(J,i) instead of a(i,J) to be consistent with matrix notation that i is along x , hence column number and J as row number.

here a(i,J) are the coefficients calculated based on type of scheme in terms of
Quote:
(F= rho*u*area and D= (mu*area)/(x(J,i+1)-x(J,i)) ).
similarly I have the equation for v.

I start with initial guess u, v, p and solve the above equations to calculate u*,v*.

Here the coefficients are calculated using u and v assumed initially.
once I have these , I solve the pressure error equation

Quote:
ap(I,J)* perror(I,J) = a(nb)*pprime(nb) + mass deficit term.
and then correct the pressure, velocity using perror terms.
The boundary conditions:
for velocity : u(:,1)=uinlet;

Quote:
at outlet du/dx=0, u( :,nx)=u( :,nx-1)
at walls u_wall=0, u(1,: )=-u(2,: ) and u(ny+1,: )=-u(ny,: )
for v velocity, v(1,: ) and v(ny,: )=0 due to wall, vin=0, which means
Quote:
v( :,1)=-v( :,2)
v, similarly at v
Quote:
at outlet =0, v(:, nx+1)=-v(:,nx )
P boundary conditions.
Quote:
at exit P = 0. P(:,nx+1)
dp/dx = at inlet .
Hence P (:,1)=P(:,2)
In the perror equation, I am solving error for Pressure nodes from 2:nx. No error at P(:,nx+1) .

Since u velocity is known at the inlet, I am assuming error to zero at P(:,1).
The problem is after few iterations it diverges like any thing. gives NaN.
I have also attached picture of the equations that are being actually solved.
Attached Images
File Type: jpg SImprimante15032514521_0001.jpg (55.3 KB, 27 views)
File Type: jpeg 14.jpeg (27.0 KB, 11 views)
alligngr8 is offline   Reply With Quote

Old   April 3, 2015, 10:25
Default
  #60
Member
 
HouKen
Join Date: Jul 2011
Posts: 67
Rep Power: 15
houkensjtu is on a distinguished road
hi everyone, I am houkensjtu, who've also been working on this problem.
Its been a while(well, 3 years...), I just was too busy to put an eye on this and I hope I could help you guys all(hopefully its not too late.)
Also my original matlab code is lost, but I will try to somehow reproduce the program, in matlab, python or C(I dont know yet).
So for those who still working on the problem, please send me your description of your problem, to my email:

houkensjtu@gmail.com

happy hacking

best regards,
houkensjtu is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
code for SIMPLE algorithm - 2D Lid driven cavity flow problem - Collocated grid h_amooie OpenFOAM Programming & Development 1 January 22, 2022 12:33
need 3D cylindrical source code for laminar flow S. D. Ding Main CFD Forum 0 July 23, 2002 03:01
What is the Better Way to Do CFD? John C. Chien Main CFD Forum 54 April 23, 2001 09:10
fluid flow fundas ram Main CFD Forum 5 June 17, 2000 22:31
Flow visualization vs. Calculated flow patterns Francisco Saldarriaga Main CFD Forum 1 August 3, 1999 00:18


All times are GMT -4. The time now is 22:14.