CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Running, Solving & CFD

CompressibleInterFoam crashing with negative temperature values

Register Blogs Community New Posts Updated Threads Search

Like Tree1Likes
  • 1 Post By student666

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   November 12, 2016, 20:12
Default CompressibleInterFoam crashing with negative temperature values
  #1
New Member
 
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 11
rzahoor is on a distinguished road
Hello Dear Foamers,

I have been trying for quite a time to solve the compressible flow of helium and water.
my boundry conditipons are as follows,
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
location "0";
object U;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 1 -1 0 0 0 0];

internalField uniform (0 0 0);

boundaryField
{
outlet-top
{
type zeroGradient;
}
outlet-front
{
type zeroGradient;
}
inlet-gas
{
type flowRateInletVelocity;
massFlowRate constant 1.1e-09;
rhoInlet 1.63;
value uniform (0 0 0);
}
inlet-liquid
{
type flowRateInletVelocity;
volumetricFlowRate constant 4.1e-12;
value uniform (0 0 0);
}
front
{
type wedge;
}
back
{
type wedge;
}
fixedWalls
{
type fixedValue;
value uniform (0 0 0);
}
}


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

object T;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [0 0 0 1 0 0 0];

internalField uniform 300;

boundaryField
{
inlet-gas
{
type fixedValue;
value uniform 300;
}

inlet-liquid
{
type fixedValue;
value uniform 300;
}

fixedWalls
{
type zeroGradient;
}

outlet-top
{
type zeroGradient;
}

outlet-front
{
type zeroGradient;
}

front
{
type wedge;
}

back
{
type wedge;
}

}

object p_rgh;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [1 -1 -2 0 0 0 0];

internalField uniform 1e5;

boundaryField
{

inlet-liquid
{
type zeroGradient;
}

inlet-gas
{
type zeroGradient;
}

fixedWalls
{
type fixedFluxPressure;
}

outlet-top
{
type uniformFixedValue;
uniformValue tableFile;
uniformValueCoeffs
{
fileName "pout";
outOfBounds clamp;
interpolationScheme linear;
};
value $internalField;
}

outlet-front
{
type uniformFixedValue;
uniformValue tableFile;
uniformValueCoeffs
{
fileName "pout";
outOfBounds clamp;
interpolationScheme linear;
};
value $internalField;
}

back
{
type wedge;
}

front
{
type wedge;
}

}

object p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions [1 -1 -2 0 0 0 0];

internalField uniform 1e5;

boundaryField
{
inlet-liquid
{
type calculated;
value $internalField;
}

inlet-gas
{
type calculated;
value $internalField;
}

fixedWalls
{
type calculated;
value $internalField;
}

outlet-top
{
type calculated;
value $internalField;
}

outlet-front
{
type calculated;
value $internalField;
}

front
{
type wedge;
}

back
{
type wedge;
}


}

I get to have the problem of case crashing with
MAXIMUM NUMBER OF ITERATIONS EXCEEDS when temperature goes to negative.
I have tried to change the pMin value and find out that it is very sensitive to the solution.
I am very much stuck.
I have opened a new tread as I have seen there are only few already exists and are quite old.
rzahoor is offline   Reply With Quote

Old   November 14, 2016, 12:07
Default
  #2
Senior Member
 
Hesam
Join Date: Feb 2015
Posts: 139
Rep Power: 11
rapierrz is on a distinguished road
Hi Rizwan,

Two suggestion:

1.Use the buoyantPressure boundary condition for pressure at walls.

2.Use the upwind scheme for div(Phi,T)
rapierrz is offline   Reply With Quote

Old   November 14, 2016, 18:14
Default
  #3
New Member
 
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 11
rzahoor is on a distinguished road
Quote:
Originally Posted by rapierrz View Post
Hi Rizwan,

Two suggestion:

1.Use the buoyantPressure boundary condition for pressure at walls.

2.Use the upwind scheme for div(Phi,T)
Hello rapierrz

I am using of3.01 which has this alternative of pressure boundary condition to buoyantPressure you suggested as "fixedFluxPressure" which I am using. and Currently I have tried with the Upwind scheme for Temperature as well.
I have tried to use the CompressiblemultiphaseInterFoam for 2 fluids to see that it effect or no.
It currently solved the temperature problem but now When I alpha.* sums increase to something 34 which is very strange even with limited schemes as well like vanLeer01.
It has been few weeks And I have used whatever possible known to me in schemes boundary conditions.
One thing that struck to me is setting the pMin in thermoPhysical properties as chaniging these values effect the run time before crashing.
I would really appreciate if you or some one can guide me through this problem.

Thanks in advance
rzahoor is offline   Reply With Quote

Old   May 9, 2017, 13:55
Default
  #4
Senior Member
 
M. C.
Join Date: May 2013
Location: Italy
Posts: 286
Blog Entries: 6
Rep Power: 16
student666 is on a distinguished road
Hi,

I'm facing same problem too, my case is a cht.
Have you found a solution?
student666 is offline   Reply With Quote

Old   May 9, 2017, 15:17
Default
  #5
New Member
 
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 11
rzahoor is on a distinguished road
Hello,

Yes i have found the solution it was just some tweeking the problem settings.
If you can explain in more detail the problem might be able to suggest something because the solution i found is problem dependent.
Regards
Rizwan
rzahoor is offline   Reply With Quote

Old   May 9, 2017, 16:52
Default
  #6
Senior Member
 
M. C.
Join Date: May 2013
Location: Italy
Posts: 286
Blog Entries: 6
Rep Power: 16
student666 is on a distinguished road
Thx for the answer, well let's begin.
My case is a conjugate heat transfer with 3 regions: two fluids and one solid.
Code:
chtMultiRegionSimpleFoam
You can see the attached picture:
  • White region is air, inlet blue arrow is inlet, outlet blue arrow is outlet.
  • Green region is air, inlet red arrow is inlet, outlet red arrow is outlet.
  • at the inlet side of the green region a scalarSemiImplicitSourceCoeffs has been applied with 5000W (absolute).
  • case is symmetryc, from the picture the front surfaces are the symmetric Patches
Here follows BC for the green region, that's where the solver blows up.
U
Code:
    airCombInlet
    {
        type            flowRateInletVelocity;
        massFlowRate    0.00261;
        extrapolateProfile yes;
        rho             rho;
        rhoInlet        1.296;
        value           uniform ( 0 0 0 );
    }
    airCombOutlet
    {
        type            zeroGradient;
    }
    simmetryComb
    {
        type            symmetry;
    }
    interCombSolido
    {
        type            noSlip;
    }
    wallComb
    {
        type            noSlip;
    }
    brac
    {
        type            noSlip;
    }
    verm
    {
        type            noSlip;
    }
    glass
    {
        type            noSlip;
    }
}
p_rgh
Code:
dimensions      [ 1 -1 -2 0 0 0 0 ];

internalField   uniform 100000;

boundaryField
{
    ".*"
    {
        type            calculated;
        value           uniform 100000;
    }
    interCombSolido
    {
        type            fixedFluxPressure;
        value           uniform 100000;
    }
    wallComb
    {
        type            fixedFluxPressure;
        value           uniform 100000;
    }
    brac    {
        type            fixedFluxPressure;
        value           uniform 100000;
    }
    verm    {
        type            fixedFluxPressure;
        value           uniform 100000;
    }
    glass
    {
        type            fixedFluxPressure;
        value           uniform 100000;
    }
    simmetryComb
    {
        type            symmetry;
    }
    airCombInlet
    {
        type            zeroGradient;
    }
    airCombOutlet
    {
        type            fixedValue;
        value           uniform 100000;
    }
}
T
Code:
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    location    "0/domain2";
    object      T;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [ 0 0 0 1 0 0 0 ];

internalField   uniform 333;

boundaryField
{
    ".*"
    {
        type            calculated;
        value           uniform 0;
    }
    airCombInlet
    {
        type            fixedValue;
        value           uniform 333;
    }
    airCombOutlet
    {
        type            inletOutlet;
        inletValue      $internalField;
        value           uniform 333;
    }
    simmetryComb
    {
        type            symmetry;
    }
    interCombSolido
    {
        type            compressible::turbulentTemperatureCoupledBaffleMixed;
        Tnbr            T;
        kappaMethod     fluidThermo;
        value           uniform 333;
    }
    wallComb
    {
        type            fixedValue;
        value           uniform 333;
    }
    glass
    {
        type            fixedValue;
        value           uniform 333;
    }
    brac    {
        type            zeroGradient;
    }
    verm    {
        type            zeroGradient;
    }
}
fvSchemes
Code:
gradSchemes
{
    default         Gauss linear;
    grad(U)         cellLimited Gauss linear 1;
}

divSchemes
{
    default         none;
    div(phi,U)      bounded Gauss upwind;
    div(phi,h)      bounded Gauss upwind;
    div(phi,e)      bounded Gauss upwind;
    div(phi,K)      bounded Gauss upwind;
    div(phi,k)      bounded Gauss upwind;
    div(phi,epsilon) bounded Gauss upwind;
    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear limited corrected 0.33;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         limited corrected 0.33;
}
fvSolution
Code:
solvers
{
    rho
    {
        solver          PCG
        preconditioner  DIC;
        tolerance       1e-7;
        relTol          0;
    }
    p_rgh
    {
        solver           GAMG;
        tolerance        1e-7;
        relTol           0.01;

        smoother         GaussSeidel;
    }

    "(U|h|k|epsilon|G|Ii)"
    {
        solver          PBiCG;
        preconditioner  DILU;
        tolerance       1e-6;
        relTol          0.1;
    }
	
	    G
    {
        $p_rgh;
        tolerance       1e-05;
        relTol          0.1;
    }
}

SIMPLE
{
    momentumPredictor yes;
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue       100000;
    rhoMin          0.2;
    rhoMax          2;
}


relaxationFactors
{
    fields
    {
        rho             1;
        p_rgh           0.7;
    }
    equations
    {
        U               0.3;
        "(h|e)"         0.3;
        k               0.3;
        epsilon         0.3;
    }
}
Code:
 Time = 111


Solving for fluid region domain0
DILUPBiCG:  Solving for Ux, Initial residual = 0.003425795, Final residual = 6.650713e-005, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.0008796441, Final residual = 1.250549e-005, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.00116966, Final residual = 1.796399e-005, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.003933221, Final residual = 0.0001691852, No Iterations 1
Min/max T:332.9558 472.7794
GAMG:  Solving for p_rgh, Initial residual = 0.003936741, Final residual = 3.912273e-005, No Iterations 4
time step continuity errors : sum local = 0.000285284, global = 3.655421e-005, cumulative = -0.1107264
Min/max rho:0.7353847 1.044312
DILUPBiCG:  Solving for epsilon, Initial residual = 0.000611032, Final residual = 1.05125e-005, No Iterations 1
DILUPBiCG:  Solving for k, Initial residual = 0.003584877, Final residual = 4.231466e-005, No Iterations 1

Solving for fluid region domain2
DILUPBiCG:  Solving for Ux, Initial residual = 0.4614772, Final residual = 0.002864163, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.3183428, Final residual = 0.003177646, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.4175119, Final residual = 0.00325495, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.07116597, Final residual = 0.0005984709, No Iterations 1
Min/max T:-2593.727 4103.124
GAMG:  Solving for p_rgh, Initial residual = 0.09838701, Final residual = 0.0005937548, No Iterations 9
time step continuity errors : sum local = 1.725629, global = 0.06717613, cumulative = -0.04355032
Min/max rho:0.2 2
DILUPBiCG:  Solving for epsilon, Initial residual = 0.435361, Final residual = 0.006099735, No Iterations 1
DILUPBiCG:  Solving for k, Initial residual = 0.2897566, Final residual = 0.003803847, No Iterations 1

Solving for solid region domain1
DICPCG:  Solving for h, Initial residual = 0.004453149, Final residual = 3.241994e-005, No Iterations 2
Min/max T:332.9824 847.455
ExecutionTime = 1175.907 s  ClockTime = 1176 s

fieldMinMax minMax write:
    min(T) = -2593.727 at location (0.04641332 0.3257848 0.06150029) on processor 1
    max(T) = 4103.124 at location (0.04349118 0.3399352 0.05850034) on processor 1
    min(mag(U)) = 0 at location (0.1318433 0.3246574 0.002) on processor 0
    max(mag(U)) = 3975.962 at location (0.04641332 0.3257848 0.06150029) on processor 1
    min(p) = -1.280926e+007 at location (0.04640306 0.3248925 0.0615002) on processor 1
    max(p) = 7.515909e+007 at location (0.04641018 0.3235 0.06149999) on processor 1

Time = 112


Solving for fluid region domain0
DILUPBiCG:  Solving for Ux, Initial residual = 0.003406134, Final residual = 6.638096e-005, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.0008713277, Final residual = 1.24063e-005, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.001160099, Final residual = 1.77978e-005, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.003906486, Final residual = 0.0001744858, No Iterations 1
Min/max T:332.9557 474.4734
GAMG:  Solving for p_rgh, Initial residual = 0.003913791, Final residual = 1.661483e-005, No Iterations 5
time step continuity errors : sum local = 0.0001209883, global = -1.862121e-005, cumulative = -0.04356895
Min/max rho:0.7327582 1.04431
DILUPBiCG:  Solving for epsilon, Initial residual = 0.0006029822, Final residual = 1.036041e-005, No Iterations 1
DILUPBiCG:  Solving for k, Initial residual = 0.00358458, Final residual = 4.245114e-005, No Iterations 1

Solving for fluid region domain2
DILUPBiCG:  Solving for Ux, Initial residual = 0.1411693, Final residual = 0.001359373, No Iterations 1
DILUPBiCG:  Solving for Uy, Initial residual = 0.1396355, Final residual = 0.001802189, No Iterations 1
DILUPBiCG:  Solving for Uz, Initial residual = 0.1626392, Final residual = 0.001579086, No Iterations 1
DILUPBiCG:  Solving for h, Initial residual = 0.03968777, Final residual = 0.000348543, No Iterations 1

job aborted:
[ranks] message
From a rough calculation, air at the green region should reach 2400K, but it gose over 5550K!!!
Then when it cools down pressure start to oscillate and when T becomes negative, it blows up...
Attached Images
File Type: png Cattura.PNG (22.6 KB, 75 views)
student666 is offline   Reply With Quote

Old   May 9, 2017, 16:56
Default
  #7
Senior Member
 
M. C.
Join Date: May 2013
Location: Italy
Posts: 286
Blog Entries: 6
Rep Power: 16
student666 is on a distinguished road
this is the output from checkMesh for green region
Code:

Create time

Create polyMesh domain2 for time = 0

Time = 0

Mesh stats
    points:           2050580
    faces:            5593869
    internal faces:   5280363
    cells:            1773483
    faces per cell:   6.131568
    boundary patches: 8
    point zones:      0
    face zones:       0
    cell zones:       4

Overall number of cells of each type:
    hexahedra:     1665623
    prisms:        3496
    wedges:        0
    pyramids:      2823
    tet wedges:    0
    tetrahedra:    1329
    polyhedra:     100212
    Breakdown of polyhedra by number of faces:
        faces   number of cells
            6   15990
            7   27867
            8   1345
            9   41192
           12   12183
           13   2
           15   1514
           18   115
           21   4

Checking topology...
    Boundary definition OK.
    Cell to face addressing OK.
    Point usage OK.
    Upper triangular ordering OK.
    Face vertices OK.
    Number of regions: 1 (OK).

Checking patch topology for multiply connected surfaces...
                   Patch    Faces   Points                  Surface topology
            airCombInlet     1930     2064  ok (non-closed singly connected)
           airCombOutlet     4526     4811  ok (non-closed singly connected)
                braciere     7452     7819  ok (non-closed singly connected)
                   glass    20480    21290  ok (non-closed singly connected)
         interCombSolido   193991   197764  ok (non-closed singly connected)
            simmetryComb    41465    42776  ok (non-closed singly connected)
             vermiculite    16810    17264  ok (non-closed singly connected)
                wallComb    26852    28166  ok (non-closed singly connected)

Checking geometry...
    Overall domain bounding box (-1.28892e-011 0.3235 0.00199986) (0.163 0.896 0.183019)
    Mesh has 3 geometric (non-empty/wedge) directions (1 1 1)
    Mesh has 3 solution (non-empty) directions (1 1 1)
    Boundary openness (4.390286e-014 2.73226e-016 -5.239765e-017) OK.
    Max cell openness = 3.663e-016 OK.
    Max aspect ratio = 14.34331 OK.
    Minimum face area = 1.858168e-009. Maximum face area = 1.639574e-005.  Face area magnitudes OK.
    Min volume = 2.982008e-013. Max volume = 7.268819e-008.  Total volume = 0.01381664.  Cell volumes OK.
    Mesh non-orthogonality Max: 75.72828 average: 8.746606
   *Number of severely non-orthogonal (> 70 degrees) faces: 303.
    Non-orthogonality check OK.
  <<Writing 303 non-orthogonal faces to set nonOrthoFaces
    Face pyramids OK.
    Max skewness = 2.409808 OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End
What is strange is that I already ran the same domain over a coarse mesh, even with more errors in the mesh, but it ran and I "succefully" (ok let's take with care) 2600K and the solution ran with no problem at all...
student666 is offline   Reply With Quote

Old   May 9, 2017, 18:24
Default
  #8
Senior Member
 
M. C.
Join Date: May 2013
Location: Italy
Posts: 286
Blog Entries: 6
Rep Power: 16
student666 is on a distinguished road
Just adding some plots about min & max for p T and U
Attached Images
File Type: png Tlimited3000up200down.png (10.9 KB, 37 views)
File Type: png pressure.png (6.6 KB, 33 views)
File Type: png U.png (12.4 KB, 33 views)
student666 is offline   Reply With Quote

Old   May 10, 2017, 10:10
Default
  #9
Senior Member
 
M. C.
Join Date: May 2013
Location: Italy
Posts: 286
Blog Entries: 6
Rep Power: 16
student666 is on a distinguished road
I think I solved my problem by setting:
Code:
    airCombInlet
    {
        type            flowRateInletVelocity;
        massFlowRate    0.00261;
        extrapolateProfile no;
        rho             rho;
        rhoInlet        1.296;
        value           uniform ( 0 0 0 );
    }
I saw it, by checking direction of the flow at the inlet side with paraview.
If it is set to yes, U vectors (magnitude) are not perpendicular to the inlet patch and, as my heat source is at the inlet side, combination of the 2 factors raises this problem...I don't think that this is the solution, but 'till now it worked for me.

Maybe this setting is useful when you have a long pipe before entering the domain and you don't know the profile....I supposed

Regards
lth likes this.
student666 is offline   Reply With Quote

Old   November 4, 2017, 14:19
Default
  #10
New Member
 
Join Date: Sep 2017
Posts: 7
Rep Power: 8
Genji is on a distinguished road
Hello Foamers
I too have run across the problem described here. After a small time (0.2 seconds) my simulation breaks down and i get negative T. Also the max courant numbers, both interface and normal become bigger than one, while the mean is 100 to 1000 times smaller. I use compressibleInterFoam to simulate a hot gas (700K) entering from the ground, in an one block mesh, and the outlet is at one of the corners at top. It is an adaptation from a case I was using with interFoam that worked. Here are my files:

blockMesh

Code:
/

convertToMeters 1;

vertices
(
	(0 0 0)
	(10 0 0)
	(10 10 0)
	(0 10 0)
	(0 0 15)
	(10 0 15)
	(10 10 15)
	(0 10 15)
	
);

blocks
(
    hex (0 1 2 3 4 5 6 7) (40 40 60) simpleGrading (1 1 1)
);

edges
(
);

boundary
(
	
	
    floor
    {
        type wall;
        faces
        (
            
			(0 3 2 1)
        );
    }
    ceiling
    {
        type wall;
        faces
        (
            
			(4 5 6 7)
        );
    }
    fixedWalls
    {
        type wall;
        faces
        (
            (0 4 7 3)
            (2 6 5 1)
            (1 5 4 0)
            (3 7 6 2)
        );
    }
);

mergePatchPairs
(
);
topoSet (for inlet and outlet)
[CODE]

actions
(
{
name jetCells;
type cellSet;
action new;
source cylinderToCell;
sourceInfo
{
p1 (5 5 -0.01);
p2 (5 5 0.26);
radius 1;
}


}

{
name f2;
type faceSet;
action new;
source cellToFace;
sourceInfo
{
set jetCells;
option all;
}
}
{
name f2;
type faceSet;
action subset;
source boxToFace;
sourceInfo
{
box (4 4 -0.1) (6 6 0.1);

}
}



{
name f1;
type faceSet;
action new;
source boxToFace;
sourceInfo
{
box (7.9 7.9 14.9) (10 10 15.1);
//boxes ((0 0 0) (1 1 1) (10 10 10)(11 11 11));
}
}
);



U

Code:
dimensions      [0 1 -1 0 0 0 0];//kg m s K mol A cd

internalField   uniform (0 0 0);//Initially the velocity is (0 0 0) m/s

boundaryField
{
    inlet
    {
    type             turbulentInlet;
	referenceField		uniform (0 0 10.5);
	fluctuationScale  		(0.1 0.1 0.5);
	value		uniform (0 0 10.5); //fixed inlet velocity
    }

    outlet
    {

        type            pressureInletOutletVelocity;
        value           uniform (0 0 0);
    }

    floor
    {
        type            noSlip;
    }

    fixedWalls
    {
        type            noSlip;
    }
	ceiling
    {
        type            noSlip;
    }
	
}
T

Code:
dimensions      [0 0 0 1 0 0 0];

internalField   uniform 300;

boundaryField
{
    inlet
    {
        type            fixedValue;
        value           uniform 700;
    }

    outlet
    {
		type            zeroGradient;
		//type            inletOutlet;
        phi             rhoPhi;
        //inletValue      $internalField;
    }

    floor
    {
        type            zeroGradient;

    }

    fixedWalls
    {
        type            zeroGradient;
    }
	ceiling
    {
        type            zeroGradient;
    }
	
}
p_rgh

Code:
 

dimensions      [1 -1 -2 0 0 0 0];

internalField   uniform 1e5;

boundaryField
{

    inlet
    {
        type            fixedFluxPressure;
		phi             rhoPhi;
        rho             rho;
        value           uniform 1e5;
    }

    outlet
    {
        type            totalPressure;
		phi             rhoPhi;
        rho             rho;
        p0              uniform 1e5;
    }

    floor
    {
        type            fixedFluxPressure;
		phi             rhoPhi;
        rho             rho;
        value           uniform 1e5;
    }

    fixedWalls
    {
        type            fixedFluxPressure;
		phi             rhoPhi;
        rho             rho;
        value           uniform 1e5;
    }
	ceiling
    {
        type            fixedFluxPressure;
		phi             rhoPhi;
        rho             rho;
        value           uniform 1e5;
    }

}
p

Code:
dimensions      [1 -1 -2 0 0 0 0];

internalField   uniform 1e5;

boundaryField
{

    inlet
    {
        type               calculated;
        value              $internalField;
    }

    outlet
    {

        type               calculated;
        value              $internalField;
    }

    floor
    {
        type               calculated;
        value              $internalField;
    }

    fixedWalls
    {
        type               calculated;
        value              $internalField;
    }
	ceiling
    {
        type               calculated;
        value              $internalField;;
    }

}
alpha.gas (we have gas and air)

Code:

dimensions      [0 0 0 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    inlet
    {
        type            fixedValue;
	value	 	uniform 1;
    }

    outlet
    {

        type            inletOutlet;
        inletValue      uniform 0;
        value           uniform 0;
    }

    floor
    {
        type            zeroGradient;
    }

    fixedWalls
    {
        type            zeroGradient;
    }
	ceiling
    {
        type            zeroGradient;
    }

}
thermophysicalProperties

Code:
phases (gas air);

pMin            pMin [1 -1 -2 0 0 0 0] 10000;

sigma           sigma [1 0 -2 0 0 0 0] 0.07;
thermophysicalProperties.air

Code:
thermoType
{
    type            heRhoThermo;
    mixture         pureMixture;
    transport       const;
    thermo          hConst;
    equationOfState perfectGas;
    specie          specie;
    energy          sensibleInternalEnergy;
}

mixture
{
    specie
    {
        molWeight   28.9;
    }
    thermodynamics
    {
        Cp          1.005;
        Hf          0;
    }
    transport
    {
        mu          1.48e-05;
        Pr          0.7;
    }
}
thermophysicalProperties.gas

Code:
thermoType
{
  
	
	type            heRhoThermo;
    mixture         pureMixture;
    transport       const;
    thermo          hConst;
    equationOfState perfectGas;
    specie          specie;
    energy          sensibleInternalEnergy;
}

mixture
{
    specie
    {
        molWeight   44.0;
    }
    equationOfState
    {
        rho         1;
    }
    thermodynamics
    {
        Cp          1.116;
        Hf          0;
    }
    transport
    {
        mu          1e-6;
        Pr          0.713;
    }
}
transportProperties

Code:
phases (gas air);

gas 
{
    transportModel  Newtonian;//Newtonian fluid assumed
    nu              [0 2 -1 0 0 0 0] 1e-06;
    rho             [1 -3 0 0 0 0 0] 0.6;
	
}

air
{
    transportModel  Newtonian;
    nu              [0 2 -1 0 0 0 0] 1.48e-05;
    rho             [1 -3 0 0 0 0 0] 1;
	
}

sigma           [1 0 -2 0 0 0 0] 0;//surface tension of gas/air
and I am calculating this as a laminar flow.

fvSolution

Code:
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //



solvers
{
	
	"alpha.gas.*"
	{
		nAlphaCorr      2;
		nAlphaSubCycles 1;
		cAlpha          1.5;
	
		MULESCorr       yes;
		nLimiterIter    3;
	
		solver          smoothSolver;
		smoother        symGaussSeidel;
		tolerance       1e-8;
		relTol          0;
	}

     cellDisplacement
    {
        solver          GAMG;
        tolerance       1e-5;
        relTol          0;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    "rho.*|pcorr.*"
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-8;
        relTol          0;
    }

    p_rgh
    {
        solver           GAMG;
        tolerance        1e-8;
        relTol           0.05;
        smoother         DICGaussSeidel;
        nPreSweeps       0;
        nPostSweeps      2;
        cacheAgglomeration on;
        agglomerator     faceAreaPair;
        nCellsInCoarsestLevel 10;
        mergeLevels      1;
    }

    p_rghFinal
    {
        $p_rgh;
        tolerance       1e-08;
        relTol          0;
    }

    "(U|k|epsilon|T)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-08;
        relTol          0.1;
    }

    "(U|k|epsilon|T)Final"
    {
        $U;
        tolerance       1e-08;
        relTol          0;
    }

}

PIMPLE
{
    momentumPredictor no;
    nOuterCorrectors 2;
    nCorrectors      4;

    nNonOrthogonalCorrectors 1;

    correctPhi              yes;

    checkMeshCourantNo      no;
    moveMeshOuterCorrectors no;

    transonic   false;
}

relaxationFactors
{
    fields
    {
    }
    equations
    {
        ".*"    1;
    }
}

// ************************************************************************* //
fvSchemes

Code:
ddtSchemes
{
    default         Euler;
}

gradSchemes
{
    default             Gauss linear;
    grad(U)             cellLimited Gauss linear 1;
}

divSchemes
{
    default         none;
    div(rhoPhi,U)       Gauss linearUpwind grad(U);
    div(phi,alpha)      Gauss vanLeer;
    div(phirb,alpha)    Gauss linear;
    div(rhoPhi,K)       Gauss linear;
    div(rhoPhi,T)       Gauss linear;
    div(phi,thermo:rho.gas) Gauss linear;
    div(phi,thermo:rho.air) Gauss linear;
    div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
    div((phi+meshPhi),p)    Gauss linear;

    div((muEff*dev2(T(grad(U))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}
controlDict

Code:
application     compressibleInterFoam;

startFrom       latestTime;

startTime       0;

stopAt          endTime;

endTime        	15;

deltaT          0.005;

writeControl    adjustableRunTime;

writeInterval   0.05;

purgeWrite      0;

writeFormat     ascii;

writePrecision  6;

writeCompression uncompressed;

timeFormat      general;

timePrecision   6;

runTimeModifiable yes;

adjustTimeStep  yes;


maxCo           1;
maxAlphaCo      1;

maxDeltaT       1;
So, where am I going wrong with this? Any help would be greatly appreciated.

Best Regards.

Last edited by Genji; November 5, 2017 at 10:55.
Genji is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
[blockMesh] blockMesh error - Negative Volume Block adoledin OpenFOAM Meshing & Mesh Conversion 2 June 22, 2016 10:44
[blockMesh] Errors during blockMesh meshing Madeleine P. Vincent OpenFOAM Meshing & Mesh Conversion 51 May 30, 2016 10:51
Ansys CFX problem: unexpected very high temperatures in premix laminar combustion faizan_habib7 CFX 4 February 1, 2016 17:00
Problem of negative temperature hnemati Main CFD Forum 0 April 1, 2014 05:25
temperature going negative behind a bluff body Sourabh Phoenics 6 December 1, 2007 18:50


All times are GMT -4. The time now is 03:59.