|
[Sponsors] |
CompressibleInterFoam crashing with negative temperature values |
|
LinkBack | Thread Tools | Search this Thread | Display Modes |
November 12, 2016, 21:12 |
CompressibleInterFoam crashing with negative temperature values
|
#1 |
New Member
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 12 |
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. |
|
November 14, 2016, 13:07 |
|
#2 |
Senior Member
Hesam
Join Date: Feb 2015
Posts: 139
Rep Power: 11 |
Hi Rizwan,
Two suggestion: 1.Use the buoyantPressure boundary condition for pressure at walls. 2.Use the upwind scheme for div(Phi,T) |
|
November 14, 2016, 19:14 |
|
#3 | |
New Member
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 12 |
Quote:
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 |
||
May 9, 2017, 16:17 |
|
#5 |
New Member
Rizwan Zahoor
Join Date: Nov 2014
Posts: 7
Rep Power: 12 |
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 |
|
May 9, 2017, 17:52 |
|
#6 |
Senior Member
|
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
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; } } 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; } } 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; } } 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; } 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 Then when it cools down pressure start to oscillate and when T becomes negative, it blows up... |
|
May 9, 2017, 17:56 |
|
#7 |
Senior Member
|
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 |
|
May 10, 2017, 11:10 |
|
#9 |
Senior Member
|
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 ); } 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 |
|
November 4, 2017, 15:19 |
|
#10 |
New Member
Join Date: Sep 2017
Posts: 7
Rep Power: 9 |
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 ( ); [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; } } 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; } } 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; } } 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;; } } 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; } } 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; 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; } } 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; } } 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 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; } } // ************************************************************************* // 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; } 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; Best Regards. Last edited by Genji; November 5, 2017 at 11:55. |
|
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
[blockMesh] blockMesh error - Negative Volume Block | adoledin | OpenFOAM Meshing & Mesh Conversion | 2 | June 22, 2016 11:44 |
[blockMesh] Errors during blockMesh meshing | Madeleine P. Vincent | OpenFOAM Meshing & Mesh Conversion | 51 | May 30, 2016 11:51 |
Ansys CFX problem: unexpected very high temperatures in premix laminar combustion | faizan_habib7 | CFX | 4 | February 1, 2016 18:00 |
Problem of negative temperature | hnemati | Main CFD Forum | 0 | April 1, 2014 06:25 |
temperature going negative behind a bluff body | Sourabh | Phoenics | 6 | December 1, 2007 19:50 |