Tutorial pq1t2_1; Part 1
Library case 621: Simulation of a Labyrinth Flow.

Preliminary remarks

This tutorial explains how to introduce parameterisation into already existing non-parameterised Q1 files.
It uses case 621 of the Input-File Library as a first simple example. This concerns steady turbulent flow, with heat transfer, through a two-dimensional labyrinthine duct. The density, specific heat, thermal conductivity and viscosity of the fluid are supposed to be constants.

The purposes of this tutorial are to show that parameterised Q1s are easier than non-parameterised ones:

We would like to remind you that an earlier tutorial concerned with flow though a labyrinth exists. It can be accessed by clicking here.

That tutorial explained how to introduce the necessary objects to the calculation domain, set their sizes and positions manually. You can examine your Q1 file if you completed all recommendations of the tutorial in question and kept the results of your work. However, you can also study the original Q1 file directly from the Input-File Library using the PHOENICS Commander environment.

Starting the tutorial

  1. We suggest that you first create a new working directory for your work. It is convenient to have a separate directory for each new project.

  2. Then open the PHOENICS Commander.

  3. Click the 'set work dir.' button in the top row, then click 'select...' and browse for the working directory you have created for this case.

  4. Click the 'Input-File Library' button in the top row.

  5. Click the 'choose by case number' button on the left, insert the case number 621 in the case number box and click on the 'Load it' button. Your actions will result in the following picture

  6. From the list of commands, choose 'See Q1'.

  7. You will see the Q1 of this case in the Commander window. Scroll through the file paying attention to what lies below the label Group 24, which, although not so marked, we shall henceforth refer to as Group 25.

    It is there where the calculation domain, the grid and all the other objects are introduced and specified.
    As you can see the sizes and positions are specified numerically, not with the help of parameters.

  8. Let us now run the case by clicking the 'Case 621' tab first and then choosing 'Run the case' command from already familiar screen.

  9. After some short time the simulation will be made and you will get a picture like this one.

  10. Click on the 'End' button, then on the 'OK' button from the Files name screen
    and your case will be opened in the VR Viewer for further inspection of the results of the simulation.

  11. For looking over the results obtained we suggest that you run the macro by clicking on the Macro button from the control panel.

  12. Click 'OK' in the MACRO Functions window and then each time follow the instructions in the Macro Message window. You might need to move the Macro Message window to prevent its obscuring the pictures like this one.

  13. Close then the VR Viewer window clicking on the top right cross of its window to return to the PHOENICS Commander. Its window could be hidden by your other opened windows. We can now examine the result file by clicking on the corresponding tab in the bottom on the window.

  14. Load the result file into the Commander editor by clicking the tab 'Result file' in the bottom row. You will see the results presented by the line-printer plot at its end.
    
       Variable    1 = P1     2 = U1     3 = W1     4 = KE     5 = EP  
         Minval=  1.000E-10 -2.974E-03  1.000E-10  5.668E-06  4.752E-08
         Maxval=  5.361E-03  1.043E-02  5.892E-02  1.781E-04  1.656E-05
         Cellav=  3.852E-03  3.696E-03  4.956E-02  1.029E-04  9.243E-06
       Variable    6 = TEM1
         Minval=  2.022E+01
         Maxval=  2.278E+01
         Cellav=  2.042E+01
     1.00 6..2.+....+..3.+3.3.+3...+....+....+....+....+....5
          .    1     3            3 3                  4 5  2
     0.90 +       3                    3  3         4  5    1
          .    3                            3  3    5    1  .
     0.80 +          1                           5  3  3 3  3
          .  3         1                       4 1  1       .
     0.70 +    2          1                    5            +
          .                 1  1       1  1 4               .
     0.60 +                       1 1       5               +
          .                               5                 .
     0.50 +       5                    4                    +
          .    5  4  5 5          4 5  5                    .
     0.40 +               5 5  5  5    2                    +
          .  5                                              .
     0.30 +  4                      2                       +
          .                                                 .
     0.20 2          2                                      +
          .                       2                         .
     0.10 +                    2                            +
          .    6  6  6 6  6                      6  6  6 6  6
     0.00 5..6.+....+....+2.6.+6..6+6..6+.6.6+.6..+....+....+
          0   .1   .2   .3   .4   .5   .6   .7   .8   .9  1.0
     the abscissa is      ISWP.  min= 1.00E+00 max= 3.81E+02
    
    Various variables at some specified monitor point, marked by numerals from 1 to 6, are plotted as a function of the iteration number, varying, as the graph indicates, from 1 to 381. The designation of variables and their limits are given above the plot.

Splitting the existing Q1 file

Before we introduce some parametrization to the basic Q1 file of the 621 case, let us create a new working directory to distinguish it from the previous one. We shall then perform all the necessary steps, one by one; and the next step will be to split the basic Q1 file into 2 parts.
  1. In the Commander window, still in your previous working directory, open the Q1 file that you have just used in the simulation. Select everything what goes immediately after
    Group 24. Dumps For Restarts
    up to the label
    STOP
    at the very bottom of the file. This is the part of Q1 where concrete data on the objects sizes and positions are stored.

  2. Cut this data into the clipboard using the hot keys combination Ctrl+X.

  3. Click on the File button and then select New from the list of commands.

    A new file will open.

  4. Inside this new file, paste the data which you have in the clipboard using the Contr+V combination. You have split the existing Q1 file into two, the first being the necessary instructions for the simulation, and the second - all numerical (for the time being) information on sizes and positions for all the objects of the labyrinth scene.

  5. Save the new file using the 'Save as' command under the 'File' button into the new directory which you have created for this case and under a new nam, say: gp25_0.

  6. Now click on the q1 tab from the bottom bar to open it. It is now time to leave there instructions about our previous actions, i.e. to give indication that all the information about the objects sizes and positions is stored in the file with the name gp25_0 in our case. You may choose any name that is meaningful for yourself.

  7. Type the following just below Group 24.

    Within familiar markers save25begin and save25end we have introduced a new PIL command incl(file_name) which loads the file named within the brackets into the instruction stack, i.e. in our case adding the required values to the list of instructions. You may read more about this command by clicking here.

  8. Having done so, save this new Q1 into the new folder with the help of already familiar Save as command.

  9. It is now time to specify your new working directory and run the case. So click on the 'Set work dir.' button in the menu bar and select the working directory in question. Answer twice 'no' to the reminding screen note to save the changed files, Q1 and new file, in the old working directory; for we shall be now working in the new one!

  10. Now run the case clicking on the 'S/E' button on the left.

  11. When the simulation is done, click on the 'END' button and open the result file inside the Commander window.

  12. If you compare both result files, you will find them to be identical. That means that the splitting of the Q1 file and its further assembling have been done correctly.

Introducing further changes

We shall here introduce our first parametrization, i.e. specify the sizes of the calculation domain via parameters.
Let us first create a new working directory for this run.
  1. In the Commander window, still in your previous working directory, open the Q1 file that you have just used in the simulation.

  2. Using the 'Save as' command copy it to your new directory.

  3. Click on 'File' and 'Open' to open the data file into the editor window.

  4. Using the 'Save as' command copy it to your new directory under any meaningful name.

  5. Switch your working directory to a newly created one with the help of the 'Set work dir.' command in the Menu bar.

  6. Click on the 'q1' button on the left to open the Q1 in the editor window.

  7. Change the name of the data file inside the 'incl' command, if you have changed it, and save the Q1 file.

  8. Click on 'File' and 'Open' to open the new data file from your new working directory in the editor window.

  9. Type the following
    xulast=2.
    yvlast=1.
    zwlast=1.
    above the line where the domain sizes are specified.

  10. Deactivate the statement with the domain sizes moving it to the third column or further and introduce instead a new one, with the domain sizes replaced with parameters: xulast, yvlast and zwlast.
    All these actions will result in the following.

  11. It is now necessary to introduce similar changes in all objects sizes and positions that depend on the domain sizes. Because this is a rather laborious and tedious work, we have done it for you. All you will need is to copy the data file named gp25_1 from the folder phoenics/d_polis/d_tuts/pq1ts/pq1t2_1 into your present working directory;
    and do not forget to introduce the corresponding changes into the incl command of the Q1 file.

  12. Now run the Q1 file from your working directory in exactly the same way as you have done before.

  13. When the simulation is finished, open the result file in the Commander editor and compare it with the previous ones: all of them should be identical.

Exploring the parametrization advantages

If you now examine the data file for this case, gp25_1, more carefully, by opening it in the Commander editor window, you will certainly find out that the sizes and positions of all objects, which make up the scene of the labyrinth case, do depend on the calculation-domain sizes, xulast, yvlast and zwlast, these being the largest values of its X-, Y- and Z- co-ordinates.

This was not evident when these positions and sizes were expressed numerically rather than via parameters.

It would convenient and reasonable to have all the objects depend on the calculation domain so that any changes introduced to the domain will result in corresponding changes of the whole scene; and we are about to find out how to achieve this.

  1. Let us first check the effect of yvlast, i.e. the domain Y-size, on the flow, which being 2-dimensional, should be unaffected by the change.

  2. Open the data file gp25_1 in the editor window if it is not still opened, clicking on File - Open - gp25_1.

  3. The default value of yvlast is 1 m. Type '2' instead meaning that we double the domain size in the Y- direction.

  4. Save the file and run the case, by clicking on the 'S/E' icon on the left.

  5. Having made the simulation, click on the 'END' button to return to the Commander editor window.

  6. Scroll to the bottom of the RESULT file where solved variables at the monitor point are plotted as a function of the iteration number. You will see that they are identical to these of the previous run, apart from the values of temperature. The results for temperature are as shown below.
    Variable    6 = TEM1
         Minval=  2.011E+01
         Maxval=  2.140E+01
         Cellav=  2.021E+01
    i.e.half the previous values.

    The reason for this is that we did not change the heat source in the L-BLOCK object. Doubling YVLAST doubled the area of the inlet and so caused twice as much air to enter; the temperature rise was therefore halved.

    To rectify this, we need to double the heat flux, i.e. in the editor to change the line

    > OBJ,    HEAT_FLUX,    0.000000E+00, 1.000000E+02
    in the file gp25_1 to
    > OBJ,    HEAT_FLUX,    0.000000E+00, 1.000000E+02*yvlast
    When you run the Satellite and solver again, you should find that the earlier temperature distribution has been recovered.

  7. It is also possible to prove the two-dimensionality of the flow in the viewer. Click the 'Run modules' button in the tool bar and then the VRV button on the left.

  8. Click 'OK' in the File names window and you will see the simulation results of this run loaded in the viewer.

  9. Click on the Contour toggle button in the tool bar to open the picture like this one.

  10. Now move probe in Y-direction paying attention to the probe value of pressure for this screen. It will remain the same. You can perform the same test for other variables, for example, velocity or temperature; or choose another position for the probe. The result of your actions will be similar: the flow is truly two-dimensional.

  11. Let us now quit the viewer by clicking on the top right cross and return to the editor to test the effect of the other domain sizes. First go to the 'View or edit files' window and click on the 'q1' icon on the left. We shall find the familiar q1 file in the editor window. However, it is not q1 that we really want to edit, but the data file gp25_1.

  12. Open the gp25_1 file clicking on 'File, Open' and selecting this very file from the 'Open' window.

  13. Increase the X-size of the domain by a factor of two, by changing the line
    xulast=2.
    to
    xulast=4.
  14. Save this change by clicking on 'File, Save File'. Then run the case.

  15. When the run is accomplished, click on the 'End' button to see the RESULT file in the editor window. Your results should be as follows.
    
       Variable    1 = P1     2 = U1     3 = W1     4 = KE     5 = EP  
         Minval=  1.000E-10  1.000E-10 -2.704E-02  5.652E-06  4.734E-08
         Maxval=  1.395E-02  1.949E-02  1.000E-10  1.087E-05  2.689E-07
         Cellav=  1.125E-02  1.773E-02 -2.253E-02  1.053E-05  2.533E-07
       Variable    6 = TEM1
         Minval=  2.001E+01
         Maxval=  2.320E+01
         Cellav=  2.017E+01
     1.00 6....+....+4.4.+4.5.+5..5+5..5+.5.5+.5.5+.5..5.5..5
          .  2 5  5  5 5  5 2  2  2 2  2  2 2  2 2  2  1    .
     0.90 +                               1 1               +
          .  5                    1 1  1                    .
     0.80 +  4              1  1                            +
          .       1  1 1  1                                 .
     0.70 +  1 1                                            +
          .                                                 .
     0.60 +                                                 +
          .                                                 .
     0.50 +                                                 +
          .                                                 .
     0.40 +                                                 +
          .                                                 .
     0.30 +  3 3                                            +
          .       3  3                                      .
     0.20 +            3  3                                 +
          .                 3  3                            .
     0.10 +                       3 3  3                    +
          .                               3 3  3 3          .
     0.00 5..6.6..6.+6.6.+6.6.+6..6+6..6+.6.6+.6.6+.6..6.6..6
          0   .1   .2   .3   .4   .5   .6   .7   .8   .9  1.0
     the abscissa is      ISWP.  min= 1.00E+00 max= 3.81E+02
     
    The numbers have been changed from their previous values; and the line-printer plot above confirms this.

  16. Let us now return to the viewer to have a closer look at the results obtained. To do so, repeat the steps g-i. You will then see a picture like this one.

    It is obvious that the X-size of the domain as well as the labyrinth itself has been noticeably changed. The labyrinth now looks longer and thinner as compared to the previous picture.
    These changes have been realized by changing only one parameter, the domain length, and all other objects being parametrically related to the domain, have been changed accordingly.

    The temperature contours will look like these.

    The velocity vectors and contours will be like these.

    Having increased the labyrinth length while retaining its height, we thereby reduced its aspect ratio zwlast/xulast; and the effect of near-wall friction became a dominant factor leading to relative increase of the recirculation zone behind the IN-BLOCK.

Exploring the parametrization advantages further

In the previous section we have exploited one advantange of parametrization namely: making the positions and sizes of all labyrinth objects depend on the domain sizes, caused any changes introduced to the latter to result in the corresponding changes to the former.
But what should we do to so as to investigate the influence of one particular object on the flow? Can parametrization be of any use here?

The answer is positive, and all we have to do is to choose an object and cut off  its connections with the domain where necessary. We are going to show all this, taking the IN-BLOCK as an example. And we shall start, as usual, with creating a new working directory for this test.

  1. Copy the data file for this case gp25_2 from the folder Phoenics/d_polis/d_tuts/pq1ts/pq1t2_1 into your new working directory. It differs from the original previous one in the following respects.

    First, new real variables are defined and set. These are:

    These varables are then inserted into the appropriate lines lower down in the data file, thus:

    While, the X- and Z-sizes of the IN-BLOCK are expressed by the corresponding parameters, its Y-size still depends on the domain as the flow remains two-dimensional.

  2. Inside the commander editor, open the previous q1 file, introduce the name of the proper data file (gp25_2) inside the incl command and save as q1 into your new working directory.

  3. Change the working directory for the new one via 'Set work dir.' button in the tool bar and run the case clicking on the 'S/E' icon on the left.

  4. Then repeat all previous actions to open the result file in the editor window.
    
       Variable    1 = P1     2 = U1     3 = W1     4 = KE     5 = EP  
         Minval=  1.000E-10 -2.974E-03  1.000E-10  5.668E-06  4.752E-08
         Maxval=  5.361E-03  1.043E-02  5.892E-02  1.781E-04  1.656E-05
         Cellav=  3.852E-03  3.696E-03  4.956E-02  1.029E-04  9.243E-06
       Variable    6 = TEM1
         Minval=  2.022E+01
         Maxval=  2.278E+01
         Cellav=  2.042E+01
     1.00 6..2.+....+..3.+3.3.+3...+....+....+....+....+....5
          .    1     3            3 3                  4 5  2
     0.90 +       3                    3  3         4  5    1
          .    3                            3  3    5    1  .
     0.80 +          1                           5  3  3 3  3
          .  3         1                       4 1  1       .
     0.70 +    2          1                    5            +
          .                 1  1       1  1 4               .
     0.60 +                       1 1       5               +
          .                               5                 .
     0.50 +       5                    4                    +
          .    5  4  5 5          4 5  5                    .
     0.40 +               5 5  5  5    2                    +
          .  5                                              .
     0.30 +  4                      2                       +
          .                                                 .
     0.20 2          2                                      +
          .                       2                         .
     0.10 +                    2                            +
          .    6  6  6 6  6                      6  6  6 6  6
     0.00 5..6.+....+....+2.6.+6..6+6..6+.6.6+.6..+....+....+
          0   .1   .2   .3   .4   .5   .6   .7   .8   .9  1.0
     the abscissa is      ISWP.  min= 1.00E+00 max= 3.81E+02
    
    If you now compare the obtained results with those for the runs made without changes in the initial data, you will find them identical. This is the proof that introduction of new parameters into the data file did not result in any changes so far.

  5. Let us check how the IN-BLOCK height, i.e. inbszz, affects the flow in the labyrinth, and for this purpose provide the initial picture of the flow velocity vectors and contours here.
    Open the viewer via the 'Run modules' windor of the commander repeating the steps g-i of the previous section. You will get picture like this one.

  6. Now return to the commander, closing the viewer and opening again the 'View and edit files' window.

  7. There click on the q1 icon on the left to open the editor window; after that open the data file gp25_2 from your working directory.

  8. The default height setting of the IN-BLOCK is inbszz=0.65. Let us increase it to 0.85, then save the change and run the case again.

  9. When the run is made, go directly to the viewer repeating your previous actions to get the picture like this one.

    You can see that the passage between the IN-BLOCK and the labyrinth wall is narrower now, which has resulted in the the corresponding increase of velocity.

  10. Let us now increase the IN-BLOCK height still more to make it block the flow completely, i.e. make inbszz=0.9.

  11. Having repeated all the previous actions we shall get the following picture in the viewer.

    The right-hand half of the picture confirms expectations: it shows that no fluid has penetrated there, and the velocity is zero everywhere.

    The left-hand half is more puzzling; for how can there be any flow there, if there is no outlet to the right?

    The explanation is that PHOENICS has been set an impossible problem: for it is still being told by the following lines in the file gp25_2:

    > OBJ,    TYPE,        INLET
    > OBJ,    PRESSURE,     0.000000E+00
    > OBJ,    VELOCITY,     5.000000E-02, 0.000000E+00, 0.000000E+00

    that there is a finite X-direction velocity of 0.05 m/s at the INLET object.

    This led to the following lines being printed near the bottom of the RESULT file:

     Nett source of R1   at patch named: OB7      (INLET   ) = 2.675250E-02
     Nett source of R1   at patch named: OB8      (OUTLET  ) =-5.449870E-11
     pos. sum=0.026752 neg. sum=-5.44987E-11
     nett sum=0.026752
    which indicates the PHOENICS could produce a solution only with a significant "nett sum" error, equal to the variable R1 which expresses the inlet mass flow rate.

    What is to be done? The simplest solution is to introduce an additional parameter, say 'invel', to set invel equal to 0.05 as the default, but then to re-set it to zero when insbzz is greater than or equal to 0.9.

    The following lines introduced to the top of gp25_2

    and the statement specifying the inlet velocityin terms of invel

     > OBJ,    VELOCITY,     invel, 0.000000E+00, 0.000000E+00
    accomplish this.

    If you re-run the case with the above settings you will see that in the result file the 'nett sum' has become negligible.

    However, if you wish to check your results in the viewer window, you might now find them to be rather strange, in the right-hand half.

    However, we suggest that you pay attention to the velocity scale bar and recognise that the largest velocity at the outlet section, being 9.24E-11 m/s, could be nothing but a result of the inevitable round-off error.

Concluding remarks

In this tutorial you have studied how You have also tested different levels of parametrization:
  1. 'global' parametrization when the positions and sizes of all the objects that make up the scene are interrelated via the calculation domain sizes; and
  2. 'local' parametrization that enables the investigation of each independent parameter effect on the flow.
Finally you made several runs and examined the results obtained, with the help of the result files and the PHOENICS viewer environment.

In two cases, unexpected results (the lower temperatures 1 and continuity 2 errors, were explained as due to aspects of the input settings which had been over-looked, to which however PHOENICS responded faithfully.

In the second part of this tutorial,, accessed by clicking here, you will learn how to:

  1. express the results of calculation in non-dimensional form;
  2. parameterize the computation-grid distribution;
  3. enable its parameters to be read from a separate file; and
  4. thereby study the influence of grid refinement on the accuracy of the results.