code_saturne User's Forum

I read the previous post on non-standard boundary conditions, but I have to ask you guys to explain it more in depth:
 
I'm trying to determine the pressure loss over a plate heat exchanger, but I cannot figure out how to set a 'constant inflow' condition for the inlet.
I would like to simulate a fixed displacement pump that 'drives' a set volumetric flow through the plate exchanger no matter how great the back-pressure is in the exchanger.
I *really* have to admit that I get lost immediately when I try to look in the usclim.f90 file
 
What I hope to do, is to run the simulation and then extract the value at the inlet and the outlet and then find the deltaP.
The simulation is so far set up fine (with 3 m/s inlet - norm) and completes without errors, but the value I get from Entree is -0.134399677E+01 (unit?) while the calculated vol. flow should be 3m/s*(12mm)^2*pi = 1,357168e-3 m3 / s = 4,8858049 m3 / h
 
How to I approach this?
 
Thanks - and congrats on the new forum :)
 
/C
 
Ps. after one (real) second the calculation goes sour and starts producing wild results - what could be the reason for this?
 

Hi Claus,

if you want to set a constant inflow, you don't need to use usclim.f90, you can do it through the GUI.

The value you get in the listing about the mass flux is in kg/s. So if your density is around 990 kg/m3, I should say everything is correct.

Differences can also arise from the discretisation of you inlet, which I understand to be a circular entrance. If you specify the exact inlet velocity (3m/s), the calculated inlet flow will be this velocity multiplied by the discretised surface, which will probably not exactly be (12mm)^2*pi.

Alternatively, you can also specify the inlet flow in the GUI, and the code will calculate the inlet velocity by dividing it by the discretised inlet surface. So you'll get an exact value for the inlet flow, but probably not strictly 3m/s for the velocity.
... but if you mesh is fine enough, it shouldn't be a problem.
 
As for the problem you get after some time of calculation, maybe you can attach the listing for us to see more clearly what happens.
 
Eventually, if you want to extract pressure values at the inlet and outlet and calculate a pressure gradient, keep in mind that the code uses P0 colocated discretisation: the values are calculated at the centers of each cells. In particular, if you use Paraview or EnSight to read the pressure values at the boundaries, note that they are only a copy made by paraview or EnSight of the value at the cell center. If you want to estimate a pressure gradient from that, you must divide the pressure difference by the distance between the centers of the cells closest to the boundaries, and not by the distance between the two boundaries.
 
Best regards
 
    --Marc--

I can be such a r*tard sometimes, it clearly states 'Mass flow' above :| I have set the density to 1000kg/m³ initially so the mass-flow is correct.
 
I have attached the 'listings' file and as you can see the derivatives and the number of iterations required starts climbing at around 15+ time step.
I've also attached some output from Paraview that graphically represents the values
 
Up till 15th time step can I trust the output enough to use the values in further calculations - and would it be correct to use the value ~0.15 bar as the pressure loss over the heat exchanger? As far as the 'plot over line' goes, it looks like it's reached a steady state.
 
Many thanks for your reply so far :)
 
/C
 
Previously Marc Sakiz wrote:

Hi Claus,
if you want to set a constant inflow, you don't need to use usclim.f90, you can do it through the GUI.

The value you get in the listing about the mass flux is in kg/s. So if your density is around 990 kg/m3, I should say everything is correct.

Differences can also arise from the discretisation of you inlet, which I understand to be a circular entrance. If you specify the exact inlet velocity (3m/s), the calculated inlet flow will be this velocity multiplied by the discretised surface, which will probably not exactly be (12mm)^2*pi.

Alternatively, you can also specify the inlet flow in the GUI, and the code will calculate the inlet velocity by dividing it by the discretised inlet surface. So you'll get an exact value for the inlet flow, but probably not strictly 3m/s for the velocity.
... but if you mesh is fine enough, it shouldn't be a problem.
 
As for the problem you get after some time of calculation, maybe you can attach the listing for us to see more clearly what happens.
 
Eventually, if you want to extract pressure values at the inlet and outlet and calculate a pressure gradient, keep in mind that the code uses P0 colocated discretisation: the values are calculated at the centers of each cells. In particular, if you use Paraview or EnSight to read the pressure values at the boundaries, note that they are only a copy made by paraview or EnSight of the value at the cell center. If you want to estimate a pressure gradient from that, you must divide the pressure difference by the distance between the centers of the cells closest to the boundaries, and not by the distance between the two boundaries.
 
Best regards
 
    --Marc--

From a quick first look at your listing, there seems to be a problem with the turbulence. Indeed, right from the beginning you have many clippings: half the cells or sometimes close to all the cells are subjected to clipping for k. The standard clipping for k in Code_Saturne consists in changing k into -k if k is negative. It is designed to correct minor non physical values in areas of low turbulence, where k is close to 0 and numerical approximations might give you a negative value close to 0. In your case, as you can see in the "INFORMATION ON CLIPPINGS" table, the value of k before clipping is quite large, as large or even larger (in absolute value) than the maximum value.
So it is clearly a sign that someting is wrong.

The first thing to do is probably to lower your time step. Indeed, your time step is 0.05s, for a velocity around 1m/s and a domain size around 0.2m, which yields a very high Courant number. You can try to lower your time step by a factor of 100 to 1000 and see what happens.
 
It has probably nothing to do with your problem, but you can add the attached files to your SRC directory. The files are a corrective patch for a potential problem with the multigrid solver for the pressure equation. It will be integrated in a corrective 2.0-rc version to be issued very soon.

Eventually, since your calculations are quite long, you can make them faster by changing the precision for the linear solvers from 1.e-8 to 1.e-5. The default value 1.e-8 originates in the early stages of development of the code, when we wanted maximum precision. Recent tests have shown that results obtained with 1.e-5 are almost identical ... but much faster. 1.e-5 is also a quite common value in other commercial codes.

You can change this value in the "Numerical Parameters"->"Equation parameters" page of the GUI ("Solver precision", to be changed for each variable).

Best regards

  --Marc--

Okay - first of all thank you for your reply - very informative.
 
I have two follow-up questions:
About the Courant–Friedrichs–Lewy condition; if I look at the current grid and determine one node with the coord. of Z=0.185 and the next node directly below it at Z=0.183428

is my CFL calc. then the following?: (0.07m/s*0.00001s)/(0.185mm−0.183428mm) ≈ 0,44529262 - is this theoretically acceptable?

The velocity of the fluid is a LOT lower this time because I devided the vol.flow with the parallel 'domains' in the exchanger. This should mean that the flow is laminar. However, I now encounter the error of

WARNING: MESH NOT ENOUGH REFINED AT THE WALL
@ ========
@ PHASE 1
@ The mesh does not seem to be enough refined at the wall
@ to be able to run a laminar simulation.
@
@ The last time step at which too large values for the
@ dimensionless distance to the wall (yplus) have been
@ observed is the time step 2
@
@ The minimum value for yplus must be lower than the
@ limit value YPLULI = 0.23810E+01

Okay - so I up the refinement of the mesh from around 300k to 1.7M cells - and nearly maxing out my 4GB of ram. STILL I get the 'mesh not enough refined at the wall'

300k cells:

 ** BOUNDARY CONDITIONS FOR SMOOTH WALLS
---------------------------------------

------------------------------------------------------------
Phase 1 Minimum Maximum
------------------------------------------------------------
Rel velocity at the wall uiptn : 0.00000E+00 0.10154E+01
Friction velocity uet : 0.13894E-02 0.10924E+00
Friction velocity uk : 0.00000E+00 0.98735E-01
Dimensionless distance yplus : 0.24844E+00 0.78270E+02
------------------------------------------------------
Nb of reversal of the velocity at the wall : 0
Nb of faces within the viscous sub-layer : 4214
Total number of wall faces : 150476
------------------------------------------------------------

1.7M cells:

** BOUNDARY CONDITIONS FOR SMOOTH WALLS
---------------------------------------

------------------------------------------------------------
Phase 1 Minimum Maximum
------------------------------------------------------------
Rel velocity at the wall uiptn : 0.00000E+00 0.11832E+01
Friction velocity uet : 0.19275E-02 0.14180E+00
Friction velocity uk : 0.00000E+00 0.12441E+00
Dimensionless distance yplus : 0.25879E+00 0.50479E+02
------------------------------------------------------
Nb of reversal of the velocity at the wall : 0
Nb of faces within the viscous sub-layer : 31294
Total number of wall faces : 583022
------------------------------------------------------------

Conversely, if I apply the K-omega SST turbulence model I get a warning of a TOO refined mesh

Now, I know this forum and the people on it are not my fluid dynamics teachers, but I hope you can spot any errors just by glancing at the numbers.

Thanks in advance.
 
/C

Hi Claus,

concerning the CFL number, your estimation is right. It is not calculated strictly like this, but roughly amounts to the same. Note that the equations are solved implicitely in Code_Saturne (at least convection and diffusion), so there is no strict CFL limitation related to stability. But since we use a predictor/corrector scheme for Navier-Stokes (SIMPLE_C scheme, velocity prediction followed by pressure correction), we do have some limitations on the time step in order to avoid velocity/pressure decoupling. We usually translate this limitation into a CFL condition. Basically we usually run our simulations with a CFL around 1 to 5 .... but it is not a strict limitation.

Moreover, if your CFL or Fourier number is too large, the linear solver (basic Jacobi for velocity) can have a hard time converging.
 
As for the wall conditions, the code uses a wall function if the y+ of the first cell is above a certain limit YPLULI (=1/kappa=2.38) and direct no-slip conditions if y+ is below YPLULI.

When you use a high Reynolds turbulence model (k-epsilon family or Rij-epsilon family), the results given if y+<YPLULI are usually not good, so the code warns you about it. It might not be a problem if only a few boundary faces are concerned, though.

When you don't use a turbulence model (laminar flow), you should not have cells where y+>YPLULI (the velocity profile is closer to parabolic than logarithmic in this area). So the code warns you about it.
 
If your problem is really laminar, it is necessary to be below YPLULI. But you don't need to refine the mesh uniformly and isotropically. You just need to refine close to the wall and in the wall-normal direction (Code_Saturne can deal with some lenght/height ratios of 10 near the wall). You can also enlarge YPLULI a bit (specifying YPLULI(1) in usini1.f90). Indeed the 1/kappa value used by Code_Saturne is not universal and another value widely used is YPLULI=10.88, which can give you some latitude and should still give you good results.

As for the k-omega SST, you shouldn't get a message, since the model is valid both close to the wall and far from the wall. However the model was designed for aerodynamics and has shown problems with heat transfer in some cases.
 
Best regards

    --Marc--

Just wanted to follow up on my case; It worked fine when I severely lowered the timestep - luckily I didn't need to track a particle through the whole geometry, otherwise the simulation would have taken 50 hours per run :S

The pressure difference I got was way too low, but that was expected since I didn't model the corrugation of the individual plates of the plate heat exchanger - otherwise the results was within the range we finally got from the supplier. The results were used to ballpark a total pressure drop of the system in the initial design phase. So mission accomplished

Attached is a screen shot from Paraview.

Thanks for the help.

/C

First post contains 1m/s inlet with 10 plates - this post 2m/s inlet with 10 plates