Nabla Logo programmer
  Advanced Search
  
  Send us a comment/query
 
  Back to main web site
 
  OpenFOAM guides
  - User guide
  - Programmer’s guide
 
  Index
 
  Changes from OpenFOAM 2.2 to 2.3
 
  Trademarks in the guides
  ©2000-2007 Nabla Ltd.

3.4 Decompression of a tank internally pressurised with water

In this example we shall investigate a problem of rapid opening of a pipe valve close to a pressurised liquid-filled tank. The prominent feature of the result in such cases is the propagation of pressure waves which must therefore be modelled as a compressible liquid.

This tutorial introduces the following OpenFOAM features for the first time:

  • Mesh refinement
  • Pressure waves in liquids

3.4.1 Problem specification

Solution domain
The domain is 2 dimensional and consists of a tank with a small outflow pipe as shown in Figure 3.9
                y


             x

        100



                50       50               240
Outlet: p =  0 bar
                            Dimensions-in-mm---
                     10    Note: image is rotated through  -90o from  normal
                                 orientation of horizontal x -axis
\special {t4ht=

Figure 3.9: Geometry of a tank with outflow pipe


Governing equations
This problem requires a model for compressibility y  \special {t4ht= in the fluid in order to be able to resolve waves propagating at a finite speed. A barotropic relationship is used to relate density r  \special {t4ht= and pressure p  \special {t4ht= are related to y  \special {t4ht=.
  • Mass continuity
    @r
---+  \~/  • (rU) = 0
@t
           \special {t4ht=
    (3.14)

  • The barotropic relationship
    @r-   -r-
@p  = K  = y
           \special {t4ht=
    (3.15)

    where K  \special {t4ht= is the bulk modulus

  • Equation 3.15 is linearised as
    r  ~~  r0 + y (p - p0)
           \special {t4ht=
    (3.16)

    where r
 0   \special {t4ht= and p
 0   \special {t4ht= are the reference density and pressure respectively such that r(p0) = r0   \special {t4ht=.

  • Momentum equation for Newtonian fluid
    @rU
-----+  \~/  • (rUU)  -  \~/  •m\ ~/ U = -  \~/ p
 @t
           \special {t4ht=
    (3.17)

Boundary conditions
 Using FoamX the following physical boundary conditions can be set:
  • outerWall is specified the wall condition;
  • axis is specified as the symmetryPlane;
  • nozzle is specified as a pressureOutlet where p = 0  \special {t4ht= bar.
  • front and back boundaries are specified as empty.
Initial conditions
U  = 0 m/s  \special {t4ht=, p = 100 bar  \special {t4ht=.
Transport properties
 
  • Dynamic viscosity of water m = 1.0 mPa  s  \special {t4ht=
Thermodynamic properties
 
  • Density of water r = 1000 kg/m3   \special {t4ht=
  • Reference pressure p =  1 bar
 0  \special {t4ht=
  • Compressibility of water y  = 4.54  10-7 s2/m2   \special {t4ht=
Solver name
sonicLiquidFoam: a compressible sonic laminar liquid flow code.
Case name
decompressionTank case located in the $OpenFOAM_TUTORIALS/sonicLiquidFoam directory.

3.4.2 Mesh Generation

The full geometry is modelled in this case; the set of vertices and blocks are given in the mesh description file below:


1  // The OpenFOAM Project // File: blockMeshDict
2  /*
3  -------------------------------------------------------------------------------
4   =========         | dictionary
5   \\      /         |
6    \\    /          | Name:   blockMeshDict
7     \\  /           | Family: FoamX configuration file
8      \\/            |
9      F ield         | OpenFOAM version: 2.3
10      O peration     | Product of Nabla Ltd.
11      A and          |
12      M anipulation  | Email: Enquiries@Nabla.co.uk
13  -------------------------------------------------------------------------------
14  */
15  // FoamX Case Dictionary.
16  
17  FoamFile
18  {
19      version         2.0;
20      format          ascii;
21  
22      root            "";
23      case            "";
24      instance        "";
25      local           "";
26  
27      class           dictionary;
28      object          blockMeshDict;
29  }
30  
31  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32  
33  convertToMeters 0.1;
34  
35  vertices
36  (
37      (0 0 -0.1)
38      (1 0 -0.1)
39      (0 0.5 -0.1)
40      (1 0.5 -0.1)
41      (1.5 0.5 -0.1)
42      (0 0.6 -0.1)
43      (1 0.6 -0.1)
44      (1.5 0.6 -0.1)
45      (0 3 -0.1)
46      (1 3 -0.1)
47      (0 0 0.1)
48      (1 0 0.1)
49      (0 0.5 0.1)
50      (1 0.5 0.1)
51      (1.5 0.5 0.1)
52      (0 0.6 0.1)
53      (1 0.6 0.1)
54      (1.5 0.6 0.1)
55      (0 3 0.1)
56      (1 3 0.1)
57  );
58  
59  blocks
60  (
61      hex (0 1 3 2 10 11 13 12) (30 20 1) simpleGrading (1 1 1)
62      hex (2 3 6 5 12 13 16 15) (30 5 1) simpleGrading (1 1 1)
63      hex (3 4 7 6 13 14 17 16) (25 5 1) simpleGrading (1 1 1)
64      hex (5 6 9 8 15 16 19 18) (30 95 1) simpleGrading (1 1 1)
65  );
66  
67  edges
68  (
69  );
70  
71  patches
72  (
73      wall outerWall
74      (
75          (0 1 11 10)
76          (1 3 13 11)
77          (3 4 14 13)
78          (7 6 16 17)
79          (6 9 19 16)
80          (9 8 18 19)
81      )
82      symmetryPlane axis
83      (
84          (0 10 12 2)
85          (2 12 15 5)
86          (5 15 18 8)
87      )
88      patch nozzle
89      (
90          (4 7 17 14)
91      )
92      empty back
93      (
94          (0 2 3 1)
95          (2 5 6 3)
96          (3 6 7 4)
97          (5 8 9 6)
98      )
99      empty front
100      (
101          (10 11 13 12)
102          (12 13 16 15)
103          (13 14 17 16)
104          (15 16 19 18)
105      )
106  );
107  
108  mergePatchPairs
109  (
110  );
111  
112  // ************************************************************************* //

In order to improve the numerical accuracy, we shall use the reference level of 1 bar  \special {t4ht= for the pressure field. Note that both the internal field level and the boundary conditions are offset by the reference level.

3.4.3 Preparing the Run

Before we commence the setup of the calculation, we need to consider the characteristic velocity of the phenomenon we are trying to capture. In the case under consideration, the fluid velocity will be very small, but the pressure wave will propagate with the speed of sound in water. The speed of sound is calculated as:

     V~ --    V~ ------------
      1-     -----1------
c =   y =    4.54  10- 7 = 1483.2m/s.
\special {t4ht=
(3.18)

For the mesh described above, the characteristic mesh size is approximately 2 m  \special {t4ht=m  \special {t4ht= (note the scaling factor of 0.1  \special {t4ht= in the blockMeshDict file). Using

       U Dt
Co  =  -----
       Dx
\special {t4ht=
(3.19)

a reasonable time step is around             -7
Dt = 5  10   s  \special {t4ht=, giving the Co  \special {t4ht= number of 0.35  \special {t4ht=, based on the speed of sound. Also, note that the reported Co  \special {t4ht= number by the code (associated with the convective velocity) will be two orders of magnitude smaller. As we are interested in the pressure wave propagation, we shall set the simulation time to 0.25 ms  \special {t4ht=. For reference, the controlDict file is quoted below.


1  // The OpenFOAM Project // File: controlDict
2  /*
3  -------------------------------------------------------------------------------
4   =========         | dictionary
5   \\      /         |
6    \\    /          | Name:   controlDict
7     \\  /           | Family: FoamX configuration file
8      \\/            |
9      F ield         | OpenFOAM version: 2.3
10      O peration     | Product of Nabla Ltd.
11      A and          |
12      M anipulation  | Email: Enquiries@Nabla.co.uk
13  -------------------------------------------------------------------------------
14  */
15  // FoamX Case Dictionary.
16  
17  FoamFile
18  {
19      version         2.0;
20      format          ascii;
21  
22      root            "";
23      case            "";
24      instance        "";
25      local           "";
26  
27      class           dictionary;
28      object          controlDict;
29  }
30  
31  // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32  
33  applicationClass sonicLiquidFoam;
34  
35  startFrom       startTime;
36  
37  startTime       0;
38  
39  stopAt          endTime;
40  
41  endTime         0.0001;
42  
43  deltaT          5e-07;
44  
45  writeControl    timeStep;
46  
47  writeInterval   20;
48  
49  cycleWrite      0;
50  
51  writeFormat     ascii;
52  
53  writePrecision  6;
54  
55  writeCompression compressed;
56  
57  timeFormat      general;
58  
59  timePrecision   6;
60  
61  runTimeModifiable yes;
62  
63  
64  // ************************************************************************* //

3.4.4 Running the case


                                                                      Pressure, p
                                                                      (bar)

                                                                           100
                                                                            90

                                                                            80
                                                                            70

                                                                            60
                                                                            50
                                                                            40

                                                                            30
                                                                            20

                                                                            10
                                                                             0




(a) At t = 50 ms         (b) At t = 100 ms        (c) At t = 150 ms
\special {t4ht=


Figure 3.10: Propagation of pressure waves


The user can run the case and view results in dxFoam. The liquid flows out through the nozzle causing a wave to move along the nozzle. As it reaches the inlet to the tank, some of the wave is transmitted into the tank and some of it is reflected. While a wave is reflected up and down the inlet pipe, the waves transmitted into the tank expand and propagate through the tank. In Figure 3.10, the pressures are shown as contours so that the wave fronts are more clearly defined than if plotted as a normal isoline plot.

If the simulation is run for a long enough time for the reflected wave to return to the pipe, we can see that negative absolute pressure is detected. The modelling permits this and has some physical basis since liquids can support tension, i.e. negative pressures. In reality, however, impurities or dissolved gases in liquids act as sites for cavitation, or vapourisation/boiling, of the liquid due to the low pressure. Therefore in practical situations, we generally do not observe pressures falling below the vapourisation pressure of the liquid; not at least for longer than it takes for the cavitation process to occur.

3.4.5 Improving the solution by refining the mesh


                                                                      Pressure, p
                                                                      (bar)

                                                                           100
                                                                            90

                                                                            80
                                                                            70

                                                                            60
                                                                            50
                                                                            40

                                                                            30
                                                                            20

                                                                            10
                                                                             0




(a) At t = 50 ms         (b) At t = 100 ms        (c) At t = 150 ms
\special {t4ht=


Figure 3.11: Propagation of pressure waves with refined mesh


Looking at the evolution of the resulting pressure field in time, we can clearly see the propagation of the pressure wave into the tank and numerous reflections from the inside walls. It is also obvious that the pressure wave is smeared over a number of cells. We shall now refine the mesh and reduce the time step to obtain a sharper front resolution. Simply edit the blockMeshDict and increase the number of cells by a factor of 4 in the x  \special {t4ht= and y  \special {t4ht= directions, i.e. block 0 becomes (120 80 1) from (30 20 1) and so on. Run blockMesh on this file. In addition, in order to maintain a Courant number below 1, the time step must be reduced accordingly to Dt =  10-7 s  \special {t4ht=. The second simulation gives considerably better resolution of the pressure waves as shown in Figure 3.11.