[heavy_user]

Hey Foamers,

i have been trying to reproduce the "simplejet" from Sandia.gov
( http://www.sandia.gov/TNF/DataArch/ProJet.html ) .

The results I produced so far are not matching the experimental data !

I would really appreciate some help here!

I tried with 2 different Geometries.

-A 2D one, like a piece of paper
-A 3D one, like a piece of cake

I have used 2 sets of boundary conditions for each Geometry
Set A : zero-gradient for outlet and sidewall.
Set B : pressureInletOutletvelocity, pressureInletoutlet, inletoutlet (the more advanced BCs)

-> which gives the 4 cases of which i have the results compared with experimental data attached.
-For experimental data I used quads/squares. For Simulation data I used crosses/daggers.
-Each axial_distance/D has a certain colour (axial_distance/D=3.775 : green, axial_distance/D=15 : yellow, axial_distance/D=30: orange, axial_distance/D=50: red; so like a traffic-light: top to bottom - green, yellow, orange, red)

So I would expect to have some daggers hitting the squares, in case that I have good results...but...

The Results are always horrible!

-The axial velocity seems to be quite close to the experimental results.
-Radial profiles are faaaaaar away.
I have massflow through the outer boundary, which should not be there. I expected this to be the cause for the bad results, but assigning phi=0 for that boundary did
not fix the problem.

I also tried different solvers (eg, simpleFoam, reactingFoam, rhoPisoFoam) without any success..

So I must be doing something essentially wrong here!??!

Since I have figuered that some other people are having those issues it might be worth to work it out in detail.

I will try to picture it in the next posts.

I have:

- Post with inlet-conditions
- Post with BC and info for 2D-case + results (using reactingFOAM)
- Post with BC and info for 3D-case + results (using reactingFOAM)
- Post with 2D geometry-details
- Post with 3D geometry-details
- Post with results for 3D-geometry with different solver

[heavy_user]

I patched the inlet and the Internal field with the velocity-profile given by sandia-data (for r/R<smaller one I used a profile from a fully developed turbulent pipe-flow, as stated in the sandia-paper).

I also used k and epsilon patches for the inlet an internal field..

For "k" i used
k(r) = 3/2 * (factor(r) *U(r)

For epsilon i used

eps(r) = U(r)^(3/2)/ (0.0052 *0.1)


I used the k-epsilon- Model for all Simulations.


I used all of them for both geometries.

Using reactingFoam i used:
tophead for C3H8 (Y=1 for r/R<1)
tophead for 02 (Y=0.21 for r/R>1)
tophead for N2 (Y=0.79 for r/R>1)
tophead for T (500K for r/R<1, 300K for r/R>1)

[heavy_user]

The 2D-Geometries flow direction is z.
I have a symmetry-plane at the z-axes.
I have 2 empty faces.


Info and BC for case "2D-A" :

Geometie: block   (symmetryPlane), Axial 2m , Radial 0,15m, D=5,5mm.
      cells (1 100 600) , Grading (1 30 40)- smallest cell axial: 0,5mm , smallest cell radial 0,26mm

Inlet:
    - vel: profile- sandia
    - eps: profile -sandia
    - k :  profile -sandia
    - T:   Tophead (500 /300)
    - C3H8,O2,N2: Tophead
    - p:   101325;

Sidewall:



    -vel :            zeroGradient;
     
    -T :             zeroGradient;

    -p :            zeroGradient;

    -k:             zeroGradient;

    -eps:           zeroGradient;
      
    -Y,C3H8;O2 :     zeroGradient;

    -N2  :          zeroGradient;


Outlet:

    -vel :            zeroGradient;
     
    -T :             zeroGradient;

    -p :            zeroGradient;

    -k:             zeroGradient;

    -eps:           zeroGradient;
      
    -Y,C3H8;O2 :     zeroGradient;

    -N2  :          zeroGradient;

                 


InternalField:

    -Profiles: U, eps, k , Tophead T, C3H8, O2,N2
               uniform  p 101325



timestep: 2,5e-6
co(first step):  Courant Number mean: 0.00744877 max: 0.342202





Mesh:

    hexahedra:     60000 
 
    Max cell openness = 2.21652e-16 OK.
    Max aspect ratio = 54.1257 OK.
    Minumum face area = 5.49269e-08. Maximum face area = 8.78116e-05.  Face area magnitudes OK.
    Min volume = 5.49269e-11. Max volume = 8.78116e-08.  Total volume = 0.00035.  Cell volumes OK.
    Mesh non-orthogonality Max: 1.40712 average: 0.526564
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.0624992 OK.

boundary conditions for case 2D-B

Geometie: Schoko-block   (symmetryPlane), Axial 2m , Radial 0,15m, D=5,5mm.
      Zellen (1 100 600) , Grading (1 30 40)- Kleinste Zelle axial: 0,5mm , kleinste Zelle radial 0,26mm

Inlet:
    - vel: profile
    - eps: profile -sandia
    - k :  profile -sandia
    - T:   Tophead
    - C3H8,O2,N2: Tophead
    - p:   zeroGradient;

Sidewall:jpg

    -vel :            type    pressureInletOutletVelocity;
                        value   uniform (0 0 0);
     
    -T :             type            inletOutlet;
                       inletValue      380;

    -p :            type            totalPressure;
                     p0              uniform 101325;
                     U               U;
                     phi             phi;
                     rho             rho;
                     psi             none;
                     gamma           1;
                     value           uniform 101325;
      
    -k:             type            inletOutlet;
                     inletValue      uniform 0.001;
             

    -eps:           type            inletOutlet;
                      inletValue      uniform 25000;
    
      
    -Y,C3H8;O2 :     type            inletOutlet;
                           inletValue      uniform 0;
    -N  :                type            inletOutlet;
                           inletValue      uniform 1;
      


Outlet:

    -vel :          type    pressureInletOutletVelocity;
                        value   uniform (0 0 0);

    -T:              zeroGradient;
       
    -p:             type            totalPressure;
                     p0              uniform 101325;
                     U               U;
                     phi             phi;
                     rho             rho;
                     psi             none;
                    gamma           1;
                    value           uniform 101325; ->JanaF

    -k:            type            inletOutlet;
                    inletValue      uniform 0.001;
            

    -eps:           type            inletOutlet;
                      inletValue      uniform 25000;
        

    -Y,C3H8;O2 :    type            inletOutlet;
                          inletValue      uniform 0;

    -N  :               type            inletOutlet;
                          inletValue      uniform 1;
                 


InternalField:

    -Profiles: U, eps, k ,  C3H8, O2,N2
                 uniform  p 101325, uniform T,





co(first step): Courant Number mean: 0.0052683 max: 0.861124
Time = 2e-06




Mesh:

    hexahedra:     60000 
 
    Max cell openness = 2.21652e-16 OK.
    Max aspect ratio = 54.1257 OK.
    Minumum face area = 5.49269e-08. Maximum face area = 8.78116e-05.  Face area magnitudes OK.
    Min volume = 5.49269e-11. Max volume = 8.78116e-08.  Total volume = 0.00035.  Cell volumes OK.
    Mesh non-orthogonality Max: 1.40712 average: 0.526564
    Non-orthogonality check OK.
    Face pyramids OK.
    Max skewness = 0.0624992 OK.
jpg

I have tried a lot of stuff in between those versions (eg. Vsidewall = Vcowflow (parallel to wall) to get rid of flow through sidewall, other solvers, fixed value conditions....) but nothing lead to better results.


I attached pictures:

Results for 2D-A (left)

[Flowtime =1 sec !]
Flux:Flux at sidewall = -0.000617497m^3/s [-37.0498 l/min]
Flux at outlet = 0.00181007m^3/s [108.604 l/min]
Flux at inlet = -0.00117736m^3/s [-70.6416 l/min]

Plot velocity-magnitude for 2D-A @ t =1sec

Pictures for 2D-B results.
A plot of velocity magnitude for 2D-B @0.37 sec (which a could not attache due to a max of 5 files) would show eddies from the shear-layer in the domain. It does NOT reach a steady state with this BCs!! I guess this is the cause for the even worse profiles (high peak in radial profiles).

Flowtime = 0.37 sec
inlet -0.00174005m^3/s
outlet 0.000597355m^3/s
sidewall 0.00112035m^3/s

[heavy_user]

The 3D-Geometries flow direction is x.
I have 2 wedge faces

Info and BC for case "3D-A" :

geo: piece of cake (wedge). Radial:0,15m ; Axial 1m ~ 200D. D=5,2mm.
      mesh (150 1 700), grading (2 1 2). -> smallests cell radial:0.55mm , smalles cell axial : 0,8mm

Inlet:
    - vel: profile
    - eps: profile -sandia
    - k :  profile -sandia
    - T:   Tophead
    - C3H8,O2,N2: Tophead
    - p:   101325;

Sidewall:

    -vel :            zeroGradient;
     
    -T :             zeroGradient;

    -p :            zeroGradient;

    -k:             zeroGradient;

    -eps:           zeroGradient;
      
    -Y,C3H8;O2 :     zeroGradient;

    -N  :           zeroGradient;


Outlet:

    -vel :            zeroGradient;
     
    -T :             zeroGradient;

    -p :            zeroGradient;

    -k:             zeroGradient;

    -eps:           zeroGradient;
      
    -Y,C3H8;O2 :     zeroGradient;

    -N  :           zeroGradient;

                 


InternalField:

    -Profiles: U, eps, k , Tophead T, C3H8, O2,N2
               uniform  p 101325



timestep: 2,5e-6
co(first step):  Courant Number mean: 0.00744877 max: 0.342202
          




Mesh:

Overall number of cells of each type:
    hexahedra:     104300            
    prisms:        700      
 
Checking geometry...
             
    Boundary openness (-2.68859e-19 -5.38961e-16 7.2511e-16) OK.                     
    Max cell openness = 2.25165e-16 OK.                                              
    Max aspect ratio = 65.5154 OK.                                                   
    Minumum face area = 2.0926e-08. Maximum face area = 2.59139e-05.  Face area magnitudes OK.
    Min volume = 2.07199e-11. Max volume = 3.57148e-08.  Total volume = 0.000980502.  Cell volumes OK.
    Mesh non-orthogonality Max: 0 average: 0                                                          
    Non-orthogonality check OK.                                                                       
    Face pyramids OK.                                                                                 
    Max skewness = 0.330796 OK.

info and bc for case 3D-B

geo: piece of cake (wedge). Radial:0,15m ; Axial 1m ~ 200D. D=5,2mm.
      mesh (150 1 700), grading (2 1 2). -> smallests cell radial:0.55mm , smalles cell axial : 0,8mm

Inlet:
    - vel: profile
    - eps: profile -sandia
    - k :  profile -sandia
    - T:   Tophead
    - C3H8,O2,N2: Tophead
    - p:   zeroGradient;

Sidewall:

    -vel :            type    pressureInletOutletVelocity;
                        value   uniform (0 0 0);
     
    -T :             type            inletOutlet;
                inletValue      380;

    -p :            type            totalPressure;
            p0              uniform 101325;
            U               U;
            phi             phi;
            rho             rho;
            psi             none;
            gamma           1;
            value           uniform 101325;
      
    -k:             type            inletOutlet;
            inletValue      uniform 0.001;
            

    -eps:           type            inletOutlet;
            inletValue      uniform 25000;
    
      
    -Y,C3H8;O2 :     type            inletOutlet;
            inletValue      uniform 0;
    -N  :           type            inletOutlet;
            inletValue      uniform 1;
      


Outlet:

    -vel :          type    pressureInletOutletVelocity;
                        value   uniform (0 0 0);

    -T:              zeroGradient;
       
    -p:             type            totalPressure;
            p0              uniform 101325;
            U               U;
            phi             phi;
            rho             rho;
            psi             none;
            gamma           1;
            value           uniform 101325;

    -k:            type            inletOutlet;
            inletValue      uniform 0.001;
            

    -eps:           type            inletOutlet;
                inletValue      uniform 25000;
        

    -Y,C3H8;O2 :    type            inletOutlet;
            inletValue      uniform 0;

    -N  :           type            inletOutlet;
            inletValue      uniform 1;
                 


InternalField:

    -Profiles: U, eps, k , Tophead T, C3H8, O2,N2
               uniform  p 101325



timestep: 2,5e-6
co(first step):  Courant Number mean: 0.00744877 max: 0.342202
          




Mesh:

Overall number of cells of each type:
    hexahedra:     104300            
    prisms:        700      
 
Checking geometry...
             
    Boundary openness (-2.68859e-19 -5.38961e-16 7.2511e-16) OK.                     
    Max cell openness = 2.25165e-16 OK.                                              
    Max aspect ratio = 65.5154 OK.                                                   
    Minumum face area = 2.0926e-08. Maximum face area = 2.59139e-05.  Face area magnitudes OK.
    Min volume = 2.07199e-11. Max volume = 3.57148e-08.  Total volume = 0.000980502.  Cell volumes OK.
    Mesh non-orthogonality Max: 0 average: 0                                                          
    Non-orthogonality check OK.                                                                       
    Face pyramids OK.                                                                                 
    Max skewness = 0.330796 OK.

I attached pictures showing results from case 3D-A

You can see that there are still waves in the flowfield which do not show up with the setting of case 3D-B.
(for the 2D geometry it is vice versa )
I have allready tried to check if they dissapear with time. But another case still had those waves at flowtimes close to 1 sec.

flowtime 0.11 sec

massflow:

inlet -0.0105862m^3/s
outlet 0.00939795m^3/s
sidewall 0.000142198m^3/s
inlet -  outlet [ 2m^3/s] =
-.00118825
Deviation from inflow through  Inlet [%]: 
11

[heavy_user]

Results for 3 D-B

 flowtime 0.474 sec

massflow:

inlet -0.00833534m^3/s
outlet 0.00740893m^3/s
sidewall 0.000926416m^3/s

I also attached a plot of velocity after 0.15 sec, which looks smooth like it should..

[heavy_user]

2D- geometry.

Mainflow is in z-Direction.
Radialdirection = y-Direction

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  1.6                                   |
|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;




vertices        
(    
    (0     0    0)  
    //(0 r*cos(2.5°) r*sin(2.5°) 0)
    (0 0.14985723323727900000 -0.00654290810480040000 )
    (0 0.14985723323727900000  0.00654290810480040000 )  
    //some point for collapsing
    (0 0  0.05)

    (1     0    0)  
    //(0 r*cos(2.5°) r*sin(2.5°) 0)
    (1 0.14985723323727900000 -0.00654290810480040000 )
    (1 0.14985723323727900000  0.00654290810480040000 )  
    //some point for collapsing
    (1 0  0.05)

);

edges
(

arc 1 2 (0 0.15 0 )

arc 5 6 (2 0.15 0)

);


blocks          
(   
    hex (0 1 2 0 4 5 6 4)      (150 1 700) simpleGrading (2 1 2)
);



patches         
(


    patch inlet
    (
        (0 2 1 0)
    )

   
    wall sidewall 
    (
       (2 6 5 1)
    )


    wedge wedge
    (
    (0 4 6 2)
    (1 5 4 0)
    )

    patch outlet    
    (
    (4 5 6 4)
    )
);

mergePatchPairs 
(
);

// ************************************************************************* //

[heavy_user]

3D- Geometry.

Mainflow is along x-Axes.
Radial flow is along y-axes.

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  1.6                                   |
|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

convertToMeters 1;




vertices        
(    
    (0     0    0)  
    //(0 r*cos(2.5°) r*sin(2.5°) 0)
    (0 0.14985723323727900000 -0.00654290810480040000 )
    (0 0.14985723323727900000  0.00654290810480040000 )  
    //some point for collapsing
    (0 0  0.05)

    (1     0    0)  
    //(0 r*cos(2.5°) r*sin(2.5°) 0)
    (1 0.14985723323727900000 -0.00654290810480040000 )
    (1 0.14985723323727900000  0.00654290810480040000 )  
    //some point for collapsing
    (1 0  0.05)

);

edges
(

arc 1 2 (0 0.15 0 )

arc 5 6 (2 0.15 0)

);


blocks          
(   
    hex (0 1 2 0 4 5 6 4)      (150 1 700) simpleGrading (2 1 2)
);



patches         
(


    patch inlet
    (
        (0 2 1 0)
    )

   
    wall sidewall 
    (
       (2 6 5 1)
    )


    wedge wedge
    (
    (0 4 6 2)
    (1 5 4 0)
    )

    patch outlet    
    (
    (4 5 6 4)
    )
);

mergePatchPairs 
(
);

// ************************************************************************* //

[heavy_user]

Using rhoPisoFoam with the settings from case 3D-B gets me the same mess.

Flowtime = 0.101 sec
(intitial unsteadynes has passed the areas of sampling)

This horrible values might be explained by the, not yet reached, steady state:

inlet -0.0083681m^3/s
outlet 0.00493538m^3/s
sidewall 0.00343236m^3/s
inlet - outlet [ 1m^3/s] =
-.00343272
deviation from flow through  Inlet [%]: 
41

[heavy_user]

Does anybody have ANY ideas what I could do to get results that have something in common with reality ??

[Vadims.geza]

Hi!

Sorry, maybe I do not understand You problem well, but why You do not use at sidewalls type fixedValue in case of no slip and type slip in case of slip? In my opinion zerogradient and pressureInletOutletVelocity are more suitable for outlet.

[heavy_user]

Originally Posted by Vadims.geza

Hi!

Sorry, maybe I do not understand You problem well, but why You do not use at sidewalls type fixedValue in case of no slip and type slip in case of slip? In my opinion zerogradient and pressureInletOutletVelocity are more suitable for outlet.

Hi Vadims.geza,

thanks for your reply!

My problems are: 1.) I have flux through sidewall, where there should
NOT be one.
2.) The results from Simulation do NOT match the
experimental data.

( The experimental data (boxes) and the simulation data (crosses) have the same colour for different hights (eg. radial_distance/D = 30 is orange ,radial_distance/D = 50 is red), so i should have red crosses in red boxed, which i have not )

I have used 2 sets of boundary conditions for each Geometry (2d-Geometry and 3D-geometry).
Set A : zero-gradient for outlet and sidewall.
Set B : pressureInletOutletvelocity, pressureInletoutlet, inletoutlet (the more advanced BCs)

So if I understand your answer correctly, I have used everything but "type slip" which should roughly the same as setting the wall-velocity=velocity_of_coflow (right?). But i will give this version a go also and get back to you with the results.

regards & nice weekend!!

[heavy_user]

Hi There,

by now I found some of my problems. I would like to share the solutions here.

First of all "slip" did a great job! THX for that advice!
Even though I still dont know why 0-gradient did not work out...

After that I figured that most essential Problem was a non constant massflow at the outlet.
(picture attached)


I tried a bunch of solvers and all of them had the non-decaying oscillation of masslfux at the outlet. The amplitude was different though.
The ONLY way to have a decaying amplitude of the massflux oscillation at the outlet I found was "wavetransmissive" as BC. This was not strait forward for me, since I thought this BC was only needed for transonic or supersonic flows...However it works.

Regards

[wayne14]

Originally Posted by heavy_user

Hi There,

by now I found some of my problems. I would like to share the solutions here.

First of all "slip" did a great job! THX for that advice!
Even though I still dont know why 0-gradient did not work out...

After that I figured that most essential Problem was a non constant massflow at the outlet.
(picture attached)


I tried a bunch of solvers and all of them had the non-decaying oscillation of masslfux at the outlet. The amplitude was different though.
The ONLY way to have a decaying amplitude of the massflux oscillation at the outlet I found was "wavetransmissive" as BC. This was not strait forward for me, since I thought this BC was only needed for transonic or supersonic flows...However it works.

Regards

Hi,
Sorry to start this old post, but I am trying to solve a jet flow too (CH4 into the air). It seems that you have solved the problem. Do you have any publications related to this?

Thank you!
Yan