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

Compressible flow solver diverging.

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 7, 2023, 03:47
Unhappy Compressible flow solver diverging.
  #1
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Recently, I have started developing my own CFD solver (hobby). I have had excellent success with in-compressible solver. (https://fluiddynamicscomputer.blogsp...-dynamics.html + https://fluiddynamicscomputer.blogsp...ierstokes.html + https://youtu.be/8RLIAkRKH64 + https://youtube.com/shorts/pKK5L6msmUI?feature=share etc.)

More recently, I have started reading about compressible equations and solvers. I am using continuity, convection diffusion (temperature), ideal gas (pressure), x and y momentum equations in time loop. I use MATLAB to test and then convert to C++ (open-source version) for publishing (speed difference is amazing, but this isn't the point of this thread ).

u(i, j) = (rho(i, j) .* u(i, j)...
- dt / (2 * h) * (rho(i+1, j) .* u(i+1, j) .* u(i+1, j) - rho(i-1, j) .* u(i-1, j) .* u(i+1, j)) ...
- dt / (2 * h) * (rho(i, j+1) .* u(i, j+1) .* v(i, j+1) - rho(i, j-1) .* u(i, j-1) .* v(i, j-1))...
- dt / (2 * h) * (p(i+1, j) - p(i-1, j)) ...
+ nu * (dt / h^2 * (u(i+1, j) - 2 * u(i, j) + u(i-1, j)) + dt / h^2 * (u(i, j+1) - 2 * u(i, j) + u(i, j-1)))) ./ rho(i, j); % x-momentum

This is supposed to be discretized x momentum in conservative form. Wonder what am I doing wrong? Tried different schemes but solution diverges. I get oscillations of velocity (horizontal lines in the domain). Doesn't matter the time step or mesh size.
fadoobaba is offline   Reply With Quote

Old   October 7, 2023, 09:24
Default
  #2
Senior Member
 
andy
Join Date: May 2009
Posts: 272
Rep Power: 18
andy_ is on a distinguished road
Quote:
Originally Posted by fadoobaba View Post
Wonder what am I doing wrong?
You haven't introduced a mechanism to control pressure-velocity decoupling at a guess. That is, there is nothing in your scheme to oppose the oscillation you see from growing to any size. There are various schemes to control it though almost all end up doing pretty much the same thing and smoothing the pressure field albeit by different amounts.
andy_ is offline   Reply With Quote

Old   October 7, 2023, 09:29
Question
  #3
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Quote:
Originally Posted by andy_ View Post
You haven't introduced a mechanism to control pressure-velocity decoupling at a guess.
Yes, no I dont think I have done pressure velocity de coupling. How to this? I will read the literature for sure, but can you guide where to start?

Currently code inside time loop has 5 equations:

Continuity for density, convection diffusion for temperature, ideal gas for pressure, u and v momentum for velocities.
fadoobaba is offline   Reply With Quote

Old   October 7, 2023, 11:32
Default
  #4
Senior Member
 
andy
Join Date: May 2009
Posts: 272
Rep Power: 18
andy_ is on a distinguished road
Quote:
Originally Posted by fadoobaba View Post
Yes, no I dont think I have done pressure velocity de coupling. How to this? I will read the literature for sure, but can you guide where to start?

Currently code inside time loop has 5 equations:

Continuity for density, convection diffusion for temperature, ideal gas for pressure, u and v momentum for velocities.
Every reasonable CFD book will discuss the range of issues that prevents a straightforward discretization of the governing transport equations for mass, momentum and energy from working. There are a range of approaches for handling the issues related to pressure-velocity decoupling, nonlineariy, compressiblity, turbulence closure and much more. Most text books will tend to focus on approaches relevant to a particular kind of flow rather than all flows. I would suggest finding a text book or course notes that is aligned with your area of interest and reading it for a few days while resisting the urge to do stuff on the computer. The links section of the website has some online sources if you don't have access to a library with a selection of CFD text books.
andy_ is offline   Reply With Quote

Old   October 7, 2023, 11:42
Default
  #5
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,782
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
I would add that a conservative FV scheme is not implemented this way! You have to compute the fluxes.
FMDenaro is offline   Reply With Quote

Old   October 7, 2023, 11:45
Thumbs up
  #6
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
I would add that a conservative FV scheme is not implemented this way! You have to compute the fluxes.
Thanks! I have read Hoffman CFD till chapter 8. I just used that knowledge to implement this. I will find another book
fadoobaba is offline   Reply With Quote

Old   October 9, 2023, 13:09
Default
  #7
New Member
 
Maksim
Join Date: Jul 2020
Posts: 10
Rep Power: 5
Eicosane is on a distinguished road
Quote:
Originally Posted by fadoobaba View Post
This is supposed to be discretized x momentum in conservative form. Wonder what am I doing wrong? Tried different schemes but solution diverges. I get oscillations of velocity (horizontal lines in the domain). Doesn't matter the time step or mesh size.
If I were you, in first place i would simplify solved system to Euler equations, i.e. exclude all diffusive terms from momentum and energy equations, and try to achieve stability and convergence for it. When i wrote my compressible solver, i used Toro's "Riemann Solvers and Numerical Methods for Fluid Dynamics". And only after that return to Navier-Stokes.

Second, as i see, your scheme is second order in space and first order in time. From my experience, this itself can lead to oscillations. Your scheme must be whether first or second order in space and time simultaneously. At least for convective term. For diffusion, i do not know, is difference in approximation order so important.

And third, looking at your equation, i have a question. Is this density: "./ rho(i, j)", at the end of equation the same as all other previous? If you solve compressible equations, it should be density from new time layer.
Eicosane is offline   Reply With Quote

Old   October 9, 2023, 14:52
Default
  #8
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Thanks for the book recommendation

In the in compressible flow part, I tried central and upwind for convective. To my surprise as well, central was more stable so I used that here. Same reason for viscous terms

Tried dividing by rho(i, j) here and rho_old(i, j) elsewhere and it showed no major difference. When I wrote in-compressible code, I used p_old (for pressure poisson), u_old, v_old etc. for velocities but then I deleted these _old matrices and replaced them with current time step one's. Code now runs 5 - 10% faster.

Here is compressible x-momentum I use:

u(i, j) = u(i, j) - u(i, j) * dt / (2 * h) .* (u(i+1, j) - u(i-1, j)) - v(i, j) * dt / (2 * h) .* (u(i, j+1) - u(i, j-1)) - dt / (2 * rho * h) * (p(i+1, j) - p(i-1, j)) + nu * (dt / h^2 * (u(i+1, j) - 2 * u(i, j) + u(i-1, j)) + dt / h^2 * (u(i, j+1) - 2 * u(i, j) + u(i, j-1))); % x-momentum

all u 's
fadoobaba is offline   Reply With Quote

Old   October 9, 2023, 15:02
Default
  #9
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,782
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by fadoobaba View Post
Thanks for the book recommendation

In the in compressible flow part, I tried central and upwind for convective. To my surprise as well, central was more stable so I used that here. Same reason for viscous terms

Tried dividing by rho(i, j) here and rho_old(i, j) elsewhere and it showed no major difference. When I wrote in-compressible code, I used p_old (for pressure poisson), u_old, v_old etc. for velocities but then I deleted these _old matrices and replaced them with current time step one's. Code now runs 5 - 10% faster.

Here is compressible x-momentum I use:

u(i, j) = u(i, j) - u(i, j) * dt / (2 * h) .* (u(i+1, j) - u(i-1, j)) - v(i, j) * dt / (2 * h) .* (u(i, j+1) - u(i, j-1)) - dt / (2 * rho * h) * (p(i+1, j) - p(i-1, j)) + nu * (dt / h^2 * (u(i+1, j) - 2 * u(i, j) + u(i-1, j)) + dt / h^2 * (u(i, j+1) - 2 * u(i, j) + u(i, j-1))); % x-momentum

all u 's

You are using a FD discretization on the quasi-linear form, I strongly discourage to use this form.

Furthermore:
- In the RHS, You are working on the same matrix u(i,j) of the update, this is wrong, you have to use two distinct arrays. What you are doing is somehow a Seidel step in time.

- in compressible flow, the diffusive term is NOT the laplacian.


I suggest to stop coding and studying the theory

- rho must be the value at (i,j)
FMDenaro is offline   Reply With Quote

Old   October 9, 2023, 15:13
Thumbs up
  #10
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Correction from previous reply:

Here is in-compressible x-momentum I use:

u(i, j) = u(i, j) - u(i, j) * dt / (2 * h) .* (u(i+1, j) - u(i-1, j)) - v(i, j) * dt / (2 * h) .* (u(i, j+1) - u(i, j-1)) - dt / (2 * rho * h) * (p(i+1, j) - p(i-1, j)) + nu * (dt / h^2 * (u(i+1, j) - 2 * u(i, j) + u(i-1, j)) + dt / h^2 * (u(i, j+1) - 2 * u(i, j) + u(i, j-1))); % x-momentum

I am learning finite volume method from Malalasekera CFD book. This blog code [https://fluiddynamicscomputer.blogsp...lab-code.html] is just for my digital CV.

Indeed the current code is finite difference as I have studied Hoffman CFD till chapter 8 completely. I will read about compressible CFD and then resume coding.
fadoobaba is offline   Reply With Quote

Old   October 9, 2023, 15:14
Default
  #11
New Member
 
Join Date: Oct 2017
Posts: 18
Rep Power: 8
fadoobaba is on a distinguished road
Quote:
Originally Posted by FMDenaro View Post
Furthermore:
- In the RHS, You are working on the same matrix u(i,j) of the update, this is wrong, you have to use two distinct arrays. What you are doing is somehow a Seidel step in time.
2 distinct arrays for each PDE makes code 5-10% slower with negligible difference in accuracy, hence I stopped using 2 arrays per PDE.
fadoobaba is offline   Reply With Quote

Old   October 9, 2023, 15:24
Default
  #12
Senior Member
 
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 6,782
Rep Power: 71
FMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura aboutFMDenaro has a spectacular aura about
Quote:
Originally Posted by fadoobaba View Post
2 distinct arrays for each PDE makes code 5-10% slower with negligible difference in accuracy, hence I stopped using 2 arrays per PDE.
And that is simply wrong for a transient computation. Your code is valid only for the steady condition.
FMDenaro is offline   Reply With Quote

Old   October 9, 2023, 18:24
Default
  #13
Senior Member
 
Join Date: Oct 2011
Posts: 242
Rep Power: 16
naffrancois is on a distinguished road
As other suggested before, CFD is not just a matter of coding finite differences disregarding the physics, the maths and stability constraints. I recommend you, for compressible flows, to consider the conservative form of the equations and a finite volume formulation. Conservative because there are strong theorems to deal with the inherent non-linearities of compressible flows (shocks, contacts) and finite volume because it is the most natural and simplest way of ensuring conservative nature of the equations. Here are two references :

- Toro: "Riemann Solvers and Numerical Methods for Fluid Dynamics" as suggested by Eicosane, for an introduction to the Riemann solution of the Euler equations and its approximations.
- Blazek: "Computational Fluid Dynamics: Principles and Applications" for a practical introduction to finite volume methods on structured and unstructured meshes
naffrancois is offline   Reply With Quote

Reply

Tags
compreesible codes, compressible, diverging, navier stokes, solver

Thread Tools Search this Thread
Search this Thread:

Advanced Search
Display Modes

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
Ncrit for a glider Xfoil. How to use it. GPT4 answer AlanMattanó Main CFD Forum 0 April 10, 2023 12:16
fluent divergence for no reason sufjanst FLUENT 2 March 23, 2016 16:08
compressible unsteady high Re and Temperature Flow: which solver to use? loreasr OpenFOAM Running, Solving & CFD 0 November 25, 2015 09:15
Solver selection for a compressible flow problem houkensjtu OpenFOAM Running, Solving & CFD 1 July 31, 2015 06:53
Solver for steady state compressible multispecies non-reacting flow? inf.vish OpenFOAM Running, Solving & CFD 0 October 1, 2013 01:09


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