[m9819348]

Hi,

I am working on some flow simulations on a micron scale (1e-6 metres).

Using the icoFoam solver on my intended model, the courant number exploded after 3 or 4 iterations. So I decided to simplify the case.

I now work on a simple rectangular pipe, with dimensions (X,Y,Z)=(0.00048, 0.00048, 0.002) metres. The mesh is a completely structured hexagonal mesh.

I have defined an inlet velocity (at Z=0.002), a pressure outlet (at Z=0.0) and normal wall boundary conditions for every other boundary. The inlet velocity is (0, 0, -0.001) m/s. The pressure in the pipe and at the outlet is 1e5 Pa (1 bar). The same problem as described above occured. The courant number explodes after only a few iterations.

If I know scale my pipe model up to metres (multiply by 1e6) and maintain the flow characteristics (same Re number, so the inlet velocity is now 1000 m/s instead of 0.001 m/s), the simulation runs without any problems!

I have tested:
- various boundary condition settings (pressure inlet, inletoutlet, ...)
- various values for the present boundary conditions (higher velocity, higher/lower pressures, ...)
- different time step sizes
- different schemes (linear, upwind)

The problem however remained the same. Below you can see a print out of the residuals of the first 10 iterations. These residuals are practically exactly the same for all the changes I tried.

So I suspect a basic error. Either in my setup (as described above) or in OpenFoam. Because the simulation (and its errors) seems to be independent of any boundary condition changes or solver setting changes.

==============================================
;Starting time loop

Time = 1e-05

Mean and max Courant Numbers = 0 0.001
BICCG: Solving for Ux, Initial residual = 0.999887, Final residual = 9.59063e-07, No Iterations 38
BICCG: Solving for Uy, Initial residual = 0.9999, Final residual = 8.2079e-07, No Iterations 38
BICCG: Solving for Uz, Initial residual = 1, Final residual = 9.24364e-07, No Iterations 39
ICCG: Solving for p, Initial residual = 1, Final residual = 8.49109e-07, No Iterations 108
time step continuity errors : sum local = 4.24544e-12, global = 1.25536e-14, cumulative = 1.25536e-14
ICCG: Solving for p, Initial residual = 0.00424226, Final residual = 8.08165e-07, No Iterations 86
time step continuity errors : sum local = 9.63602e-10, global = 1.59232e-12, cumulative = 1.60488e-12
ExecutionTime = 12.92 s ClockTime = 13 s

Time = 2e-05

Mean and max Courant Numbers = 0.000171858 0.00108216
BICCG: Solving for Ux, Initial residual = 0.54666, Final residual = 7.87019e-07, No Iterations 34
BICCG: Solving for Uy, Initial residual = 0.54666, Final residual = 7.87019e-07, No Iterations 34
BICCG: Solving for Uz, Initial residual = 0.414357, Final residual = 9.10675e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.286481, Final residual = 8.46628e-07, No Iterations 98
time step continuity errors : sum local = 6.00845e-11, global = -2.96465e-13, cumulative = 1.30841e-12
ICCG: Solving for p, Initial residual = 0.00431205, Final residual = 9.16835e-07, No Iterations 84
time step continuity errors : sum local = 4.12726e-09, global = -2.86835e-13, cumulative = 1.02158e-12
ExecutionTime = 24.36 s ClockTime = 24 s

Time = 3e-05

Mean and max Courant Numbers = 0.000388712 0.0058744
BICCG: Solving for Ux, Initial residual = 0.793024, Final residual = 6.6019e-07, No Iterations 35
BICCG: Solving for Uy, Initial residual = 0.793024, Final residual = 6.6019e-07, No Iterations 35
BICCG: Solving for Uz, Initial residual = 0.504602, Final residual = 7.54358e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.353671, Final residual = 9.81328e-07, No Iterations 97
time step continuity errors : sum local = 2.14818e-10, global = 5.53171e-13, cumulative = 1.57475e-12
ICCG: Solving for p, Initial residual = 0.00437745, Final residual = 8.03115e-07, No Iterations 84
time step continuity errors : sum local = 1.34631e-08, global = 6.95756e-11, cumulative = 7.11503e-11
ExecutionTime = 35.24 s ClockTime = 35 s

Time = 4e-05

Mean and max Courant Numbers = 0.000984455 0.0125623
BICCG: Solving for Ux, Initial residual = 0.949905, Final residual = 7.04181e-07, No Iterations 35
BICCG: Solving for Uy, Initial residual = 0.949905, Final residual = 7.04181e-07, No Iterations 35
BICCG: Solving for Uz, Initial residual = 0.791293, Final residual = 9.88264e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.339371, Final residual = 6.88011e-07, No Iterations 97
time step continuity errors : sum local = 5.921e-10, global = -4.02724e-14, cumulative = 7.111e-11
ICCG: Solving for p, Initial residual = 0.00444162, Final residual = 8.52582e-07, No Iterations 83
time step continuity errors : sum local = 5.30509e-08, global = 3.33387e-10, cumulative = 4.04497e-10
ExecutionTime = 46.2 s ClockTime = 46 s

Time = 5e-05

Mean and max Courant Numbers = 0.0039436 0.0561303
BICCG: Solving for Ux, Initial residual = 0.924943, Final residual = 6.20131e-07, No Iterations 35
BICCG: Solving for Uy, Initial residual = 0.924943, Final residual = 6.20131e-07, No Iterations 35
BICCG: Solving for Uz, Initial residual = 0.689644, Final residual = 7.95788e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.347068, Final residual = 8.67243e-07, No Iterations 96
time step continuity errors : sum local = 2.74462e-09, global = -3.8124e-12, cumulative = 4.00684e-10
ICCG: Solving for p, Initial residual = 0.00449369, Final residual = 9.90022e-07, No Iterations 82
time step continuity errors : sum local = 2.28864e-07, global = -1.45197e-09, cumulative = -1.05128e-09
ExecutionTime = 57.27 s ClockTime = 58 s

Time = 6e-05

Mean and max Courant Numbers = 0.0141163 0.20115
BICCG: Solving for Ux, Initial residual = 0.940933, Final residual = 5.78777e-07, No Iterations 35
BICCG: Solving for Uy, Initial residual = 0.940933, Final residual = 5.78777e-07, No Iterations 35
BICCG: Solving for Uz, Initial residual = 0.716215, Final residual = 7.97297e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.348306, Final residual = 9.95881e-07, No Iterations 95
time step continuity errors : sum local = 1.17936e-08, global = 4.22655e-11, cumulative = -1.00902e-09
ICCG: Solving for p, Initial residual = 0.0045422, Final residual = 8.61181e-07, No Iterations 82
time step continuity errors : sum local = 7.39444e-07, global = 4.11725e-09, cumulative = 3.10823e-09
ExecutionTime = 67.45 s ClockTime = 68 s

Time = 7e-05

Mean and max Courant Numbers = 0.0523363 0.765233
BICCG: Solving for Ux, Initial residual = 0.945119, Final residual = 9.80385e-07, No Iterations 34
BICCG: Solving for Uy, Initial residual = 0.945119, Final residual = 9.80385e-07, No Iterations 34
BICCG: Solving for Uz, Initial residual = 0.708878, Final residual = 7.72249e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.350959, Final residual = 8.54684e-07, No Iterations 95
time step continuity errors : sum local = 3.76856e-08, global = -8.42649e-11, cumulative = 3.02397e-09
ICCG: Solving for p, Initial residual = 0.00458737, Final residual = 7.71143e-07, No Iterations 82
time step continuity errors : sum local = 2.45929e-06, global = -1.10461e-08, cumulative = -8.02209e-09
ExecutionTime = 78.99 s ClockTime = 80 s

Time = 8e-05

Mean and max Courant Numbers = 0.192457 2.86185
BICCG: Solving for Ux, Initial residual = 0.951022, Final residual = 9.26774e-07, No Iterations 34
BICCG: Solving for Uy, Initial residual = 0.951022, Final residual = 9.26774e-07, No Iterations 34
BICCG: Solving for Uz, Initial residual = 0.710894, Final residual = 7.52262e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.353132, Final residual = 7.85644e-07, No Iterations 95
time step continuity errors : sum local = 1.29109e-07, global = 1.67507e-10, cumulative = -7.85458e-09
ICCG: Solving for p, Initial residual = 0.00462957, Final residual = 9.49266e-07, No Iterations 81
time step continuity errors : sum local = 1.12486e-05, global = -3.66109e-08, cumulative = -4.44655e-08
ExecutionTime = 89.02 s ClockTime = 90 s

Time = 9e-05

Mean and max Courant Numbers = 0.710661 10.7392
BICCG: Solving for Ux, Initial residual = 0.954611, Final residual = 9.07716e-07, No Iterations 34
BICCG: Solving for Uy, Initial residual = 0.954611, Final residual = 9.07716e-07, No Iterations 34
BICCG: Solving for Uz, Initial residual = 0.710387, Final residual = 8.59814e-07, No Iterations 35
ICCG: Solving for p, Initial residual = 0.355131, Final residual = 9.05752e-07, No Iterations 93
time step continuity errors : sum local = 5.54213e-07, global = -8.5073e-09, cumulative = -5.29728e-08
ICCG: Solving for p, Initial residual = 0.00467167, Final residual = 8.40575e-07, No Iterations 81
time step continuity errors : sum local = 3.6961e-05, global = 3.36836e-07, cumulative = 2.83863e-07
ExecutionTime = 100.2 s ClockTime = 101 s

Time = 0.0001

Mean and max Courant Numbers = 2.61587 40.3566
BICCG: Solving for Ux, Initial residual = 0.958569, Final residual = 9.10853e-07, No Iterations 34
BICCG: Solving for Uy, Initial residual = 0.958569, Final residual = 9.10853e-07, No Iterations 34
BICCG: Solving for Uz, Initial residual = 0.710484, Final residual = 8.79641e-07, No Iterations 37
ICCG: Solving for p, Initial residual = 0.358534, Final residual = 9.95292e-07, No Iterations 94
time step continuity errors : sum local = 2.26458e-06, global = 7.7304e-09, cumulative = 2.91593e-07
ICCG: Solving for p, Initial residual = 0.00469782, Final residual = 8.45221e-07, No Iterations 81
time step continuity errors : sum local = 0.00013883, global = -1.07493e-06, cumulative = -7.83336e-07
ExecutionTime = 110.85 s ClockTime = 112 s
============================================

My questions are:
1. Am I making a mistake here? Am I overlooking something?
2. Does OpenFoam has a basic problem with solving meshes on a micron scale?
3. How can I solve this situation? How can I make OpenFoam run simulations on a micron scale?

Thanks,

Ries Bouwman

[hjasak]

A few comments:
- if you are running an incompressible flow solver, the pressure level at the outlet is irrelevant and you can put it to zero rather than 1 bar. It should make no difference (only round-off in solver operations), but please bear with me
- more importantly, from the log file I can see that the momentum solver is working extremely hard, which is unusual "under normal circumstances". I would have expected it to do 1-2 BiCCG sweeps. This would indicate that your problem is totally dominated by the diffusion in the momentum equation, which is again quite unusual. Please check your Reynolds number (this may be what you want).
- there should be no fundamental problems on running the code on this scale (in fact, this is still quite big) :-)

For such cases, you will be governed by the diffusion Courant number (see Patankar) rather than the convection one reported by the code. If this is so, you can try dropping the time step by 1-2 orders of magnitude and see if it makes things better.

Two final checks: you are running double precision, right? (please say you do) + writing out the mesh in binary to keep the precision of point locations. Secondly, can you tell me something about how the code actually fails, e.g. does the solution go stupid next to a boundary, slowly or rapidly, in the pressure or velocity etc.

Hrv

[m9819348]

Dear Mr Jasak,

thank you for your quick reply!

Regarding your questions:
1. Re ~ 0.5, my material parameters are those of normal water. I see no real reason why the problem should become diffusion dominated.
2. However, reducing the time step from 1e-5 to 1e-7 made the simulation stable. It now only needs about 3 BiCCG sweeps.
3. As far as I know, the version OF1.2 runs automatically in double precision (so does the 1.3). So, yes, I am running the simulation in DP.
4. I have written out my mesh in Gambit (.msh). I can not choose an option to write it binary. Can I change my mesh to a binary mesh in OF1.2?
5. Finally, the solution goes bezerk next to the wall boundary conditions, especially near the outlet. This happens rapidly (within one time step and after only 3 or 4 iterations). Both the pressure field and the velocity field show idiotic results.

I have a few more questions:
a. Where exactly can I change the precision in OF1.2? I did not find this in the manual, nor in the discussion groups.
b. Is there a more precise way of exporting a gambit mesh and then importing (converting) this mesh to OpenFoam?

[hjasak]

Hi,

Good: progress! :-)

a) In order to run single precision, you would need to re-compile the code yourself or download a single precision version. This is controlled by:

setenv WM_PRECISION_OPTION DP (or SP)

in

~/.OpenFOAM-1.3/cshrc

b) If the mesh is generated in OpenFOAM, you can run blockMesh with binary output format or change the writePrecision in system/controlDict (there's 15 significant places in DP). Unfortunately, I don't read Fluent/Gambit binary format in the converter so the limitation comes from the accuracy you can get in vertex positions in a Gambit/Fluent file... and I don't know how to control those. Once the mesh is converted, it does not move so there's no point in switching the format.

Enjoy,

Hrv

[m9819348]

Hi,

thank you for your (very!) quick response.

My run is now proceeding in a proper manner. It will not be long before I can take a look at the solution. Thank you for your help.

Unfortunately, I can not create a blockmesh for the model that actually interests me, since it is to complicated for that. But I believe that the help you offered above is enough for me to let the model run properly.

I will try that single precision solver as soon as my mesh becomes very unstructured. I believe that single precision works better then than double, doesn't it?

Thank you once more and goob-bye.

Ries

[hjasak]

I would recommend you stick to double precision whenever possible - it eliminates a lot of potential problem sources. Unless you are planning to run hundred million cell runs and relatively simple physics, there's no need to go to single precision because there's plenty of memory available.

Hrv

[eugene]

Just a thought:

1. Make the mesh in fluent in micrometers.
2. Convert it to OF format.
3. Change the write format to binary or increase the write precision to 12.
4. Use transformPoints to scale back to meters.

Though this will only help with a limited subset of problems.

[hoogland]

Hi!

Have you got a more specific Patankar reference for that diffusive Co number? I can't find an explicit discussion in "Numerical Heat Transfer and Fluid Flow". I've got a case of a high p reservoir feeding supersonically into a quiescent chamber, with the same symptoms (mom solver working hard, Co rises etc) and suspect this explanation but want to learn more, i.e. want to maximise dt. Any suggestions on best quantifying/displaying diffusive flux v convective flux contributions in OF?

Hoogs