# Problems with Turbulence Modeling

 Register Blogs Members List Search Today's Posts Mark Forums Read

 November 26, 2009, 05:48 Problems with Turbulence Modeling #1 New Member   Jignesh Join Date: Sep 2009 Posts: 6 Rep Power: 9 Hi, I need to simulate flow over a cylinder for the following cases: 1) Inviscid Flow 2) Viscous Subsonic 3) Viscous Subsonic with RAS Turbulence 4) Supersonic I'm able to simulate the first two (using potentialFoam and icoFoam respectively), but am getting this error while trying to simulate the last two (using pisoFoam and sonicFoam respectiely) Create time Create mesh for time = 0 Reading field p Reading field U Reading/calculating face flux field phi Selecting incompressible transport model Newtonian Selecting turbulence model type RASModel Selecting RAS turbulence model kEpsilon kEpsilonCoeffs { Cmu 0.09; C1 1.44; C2 1.92; sigmaEps 1.3; } Starting time loop Time = 0.0005 Courant Number mean: 0.000271663 max: 0.00134942 DILUPBiCG: Solving for Ux, Initial residual = 1, Final residual = 5.73156e-08, No Iterations 1 DILUPBiCG: Solving for Uy, Initial residual = 1, Final residual = 6.21305e-08, No Iterations 1 DICPCG: Solving for p, Initial residual = 1, Final residual = 0.0913066, No Iterations 31 time step continuity errors : sum local = 8.14153e-05, global = -4.30474e-06, cumulative = -4.30474e-06 DICPCG: Solving for p, Initial residual = 0.0302236, Final residual = 6.3797e-07, No Iterations 49 time step continuity errors : sum local = 7.42487e-08, global = 1.87645e-10, cumulative = -4.30455e-06 #0 Foam::error:rintStack(Foam::Ostream&) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so" #1 Foam::sigFpe::sigFpeHandler(int) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so" #2 Uninterpreted: #3 Foam::divide(Foam::Field&, Foam::UList const&, Foam::UList const&) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so" #4 void Foam::divide(Foam::GeometricField&, Foam::GeometricField const&, Foam::GeometricField const&) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so" #5 Foam::tmp > Foam:perator/(Foam::tmp > const&, Foam::GeometricField const&) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so" #6 Foam::incompressible::RASModels::kEpsilon::correct () in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libincompressibleRASModels.so" #7 main in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/applications/bin/linuxGccDPOpt/pisoFoam" #8 __libc_start_main in "/lib/tls/i686/cmov/libc.so.6" #9 _start at /usr/src/packages/BUILD/glibc-2.9/csu/../sysdeps/i386/elf/start.S:122 Floating point exception Note that the geometry remains the same in all the cases, only the solver and the relevant inputs and boundary types change depending upon the case. Also, for the supersonic case, the solution blows up (as shown above) after several iteration, whereas for the subsonic turbulent case, it blows up right from the start. Can anyone please gimme a clue as to what is it that I'm doing wrong? Regards, Jignesh

 November 26, 2009, 08:49 #2 Member   Ulrich Heck Join Date: Mar 2009 Location: Krefeld, Germany Posts: 41 Rep Power: 9 Hi Jignesh, I think there is something wrong with your initial turbulence conditions. This might be caused by your turbulence div-schemes. Waht do you use? Ulrich

 November 26, 2009, 10:46 #4 Member   Ulrich Heck Join Date: Mar 2009 Location: Krefeld, Germany Posts: 41 Rep Power: 9 Hi Jignesh I think you should not force k to zero at your patch "right". Maybe you should try zero gradient instead. If this doesn't help. Consider more stable schemes for turbulence, e.g. div(phi,k) Gauss upwind; div(phi,epsilon) Gauss upwind; and switch to more accurate schemes later. But in most cases, when the solution crashes so soon, there is something wrong with the initial conditions or bcs. You can also set turbulence Intensities and mixingLenght for inlets e.g. k: type turbulentIntensityKineticEnergyInlet; intensity 0.01; value \$internalField; e.g. omega: type turbulentMixingLengthFrequencyInlet; mixingLength 0.04; k k; value \$internalField; For epsilon is turbulentMixingLengthDissipationInlet I think. Hope this helps Ulli

 November 26, 2009, 16:12 #5 New Member   Jignesh Join Date: Sep 2009 Posts: 6 Rep Power: 9 Hey Ulrich, Thanks for the suggestion; it worked! i got the simulation running flawlessly by NOT forcing k to be equal to zero! now I jst need to ensure that the results that I have are indeed correct... however, the problem with the supersonic case seems to be something different, since it is blowing up suddenly after 0.0027secs: Time = 0.00273952 Courant Number mean: 0.00110348 max: 0.0129647 diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0 DILUPBiCG: Solving for Ux, Initial residual = 0.970821, Final residual = 4.52189e-06, No Iterations 27 DILUPBiCG: Solving for Uy, Initial residual = 0.950914, Final residual = 8.50896e-06, No Iterations 26 DILUPBiCG: Solving for e, Initial residual = 0.956971, Final residual = 5.481e-06, No Iterations 27 DILUPBiCG: Solving for p, Initial residual = 0.100512, Final residual = 7.64695e-09, No Iterations 13 DILUPBiCG: Solving for p, Initial residual = 0.0364326, Final residual = 9.77965e-09, No Iterations 7 DILUPBiCG: Solving for p, Initial residual = 0.00127526, Final residual = 2.52387e-09, No Iterations 9 diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0 time step continuity errors : sum local = 2.60677e-09, global = 1.43791e-09, cumulative = 2.09749e-08 DILUPBiCG: Solving for p, Initial residual = 0.0276749, Final residual = 2.05278e-10, No Iterations 15 DILUPBiCG: Solving for p, Initial residual = 0.019402, Final residual = 1.06515e-09, No Iterations 7 DILUPBiCG: Solving for p, Initial residual = 0.00143388, Final residual = 7.04595e-09, No Iterations 5 diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0 time step continuity errors : sum local = 8.05887e-08, global = -4.09812e-08, cumulative = -2.00063e-08 DILUPBiCG: Solving for epsilon, Initial residual = 0.415028, Final residual = 6.03136e-09, No Iterations 35 bounding epsilon, min: -4.23037e+26 max: 6.89055e+26 average: 2.43139e+23 DILUPBiCG: Solving for k, Initial residual = 0.980292, Final residual = 2.48219e-09, No Iterations 34 bounding k, min: -3.11137e+16 max: 2.40184e+17 average: 7.80979e+15 ExecutionTime = 6421.04 s ClockTime = 6739 s Time = 0.00273956 Courant Number mean: 0.108861 max: 259.099 diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0 DILUPBiCG: Solving for Ux, Initial residual = 0.949675, Final residual = 8.48593e-06, No Iterations 26 DILUPBiCG: Solving for Uy, Initial residual = 0.902717, Final residual = 8.30647e-06, No Iterations 26 DILUPBiCG: Solving for e, Initial residual = 0.954222, Final residual = 6.66186e-06, No Iterations 26 Maximum number of iterations exceeded#0 Foam::error:rintStack(Foam::Ostream&) in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so" #1 Foam::error::abort() in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libOpenFOAM.so" #2 Foam::ePsiThermo > > > >::calculate() in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libbasicThermophysicalModels.so" #3 Foam::ePsiThermo > > > >::correct() in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/lib/linuxGccDPOpt/libbasicThermophysicalModels.so" #4 main in "/home/jignesh/OpenFOAM/OpenFOAM-1.6/applications/bin/linuxGccDPOpt/sonicFoam" #5 __libc_start_main in "/lib/tls/i686/cmov/libc.so.6" #6 _start at /usr/src/packages/BUILD/glibc-2.9/csu/../sysdeps/i386/elf/start.S:122 From function specieThermo::T(scalar f, scalar T0, scalar (specieThermo::*F)(const scalar) const, scalar (specieThermo::*dFdT)(const scalar) const) const in file /home/dm2/henry/OpenFOAM/OpenFOAM-1.6/src/thermophysicalModels/specie/lnInclude/specieThermoI.H at line 68. FOAM aborting Aborted The Courant No is SUDDENLY shooting up by a factor of 20,000 in jst one iteration. Infact the results till t=0.0027 seem to be fairly decent when viewed in paraFoam. but I want to run the simulation till 0.01 to ensure it has reached steady-state; and I have no clue as to WHAT is causing this blowup. Any ideas about trouble shooting? Has it got anything to do with the thermophysical model? i'm reproducing the files for the supersonic case to make things more clear: Again, many thanks for the help! Regards, Jignesh ======================== FoamFile { version 2.0; format ascii; class volScalarField; object p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [1 -1 -2 0 0 0 0]; internalField uniform 100000; boundaryField { inlet { type fixedValue; value uniform 100000; } outlet { type waveTransmissive; field p; phi phi; rho rho; psi psi; gamma 1.3; fieldInf 100000; lInf 1; value uniform 100000; } topwall { type zeroGradient; } symplane { type symmetryPlane; } cylinder { type zeroGradient; } frontAndBack { type empty; } } ================================ FoamFile { version 2.0; format ascii; class volVectorField; location "0"; object U; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 1 -1 0 0 0 0]; internalField uniform (650 0 0); boundaryField { topwall { type supersonicFreestream; pInf 100000; TInf 300; UInf (650 0 0); gamma 1.28418; value uniform (650 0 0); } symplane { type symmetryPlane; } frontAndBack { type empty; } cylinder { type fixedValue; value uniform (0 0 0); } inlet { type fixedValue; value uniform (650 0 0); } outlet { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } } ================================ FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object T; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 0 0 1 0 0 0]; internalField uniform 300; boundaryField { inlet { type fixedValue; value uniform 300; } outlet { type inletOutlet; inletValue uniform 300; value uniform 300; } symplane { type symmetryPlane; } topwall { type inletOutlet; inletValue uniform 300; value uniform 300; } cylinder { type zeroGradient; } frontAndBack { type empty; } } =============================== FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object k; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -2 0 0 0 0]; internalField uniform 1000; boundaryField { inlet { type fixedValue; value uniform 1000; } outlet { type inletOutlet; inletValue uniform 1000; value uniform 1000; } topwall { type inletOutlet; inletValue uniform 1000; value uniform 1000; } symplane { type symmetryPlane; } cylinder { type compressible::kqRWallFunction; value uniform 1000; } frontAndBack { type empty; } } =================================== FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object mut; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [1 -1 -1 0 0 0 0]; internalField uniform 0; boundaryField { inlet { type calculated; value uniform 0; } outlet { type calculated; value uniform 0; } topwall { type calculated; value uniform 0; } symplane { type symmetryPlane; } cylinder { type mutWallFunction; Cmu 0.09; kappa 0.41; E 9.8; value uniform 0; } frontAndBack { type empty; } } ======================================== FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object epsilon; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -3 0 0 0 0]; internalField uniform 25000; boundaryField { inlet { type fixedValue; value uniform 25000; } outlet { type inletOutlet; inletValue uniform 25000; value uniform 25000; } topwall { type inletOutlet; inletValue uniform 25000; value uniform 25000; } symplane { type symmetryPlane; } cylinder { type compressible::epsilonWallFunction; Cmu 0.09; kappa 0.41; E 9.8; value uniform 25000; } frontAndBack { type empty; } } ================================================= FoamFile { version 2.0; format ascii; class dictionary; location "system"; object fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // ddtSchemes { default Euler; } gradSchemes { default Gauss linear; } divSchemes { default none; div(phi,U) Gauss limitedLinearV 1; div(phi,k) Gauss upwind; div(phi,epsilon) Gauss upwind; div(phi,R) Gauss upwind; div(R) Gauss linear; div(phid,p) Gauss limitedLinear 1; div(phiU,p) Gauss limitedLinear 1; div(phi,e) Gauss limitedLinear 1; div((muEff*dev2(grad(U).T()))) Gauss linear; } laplacianSchemes { default none; laplacian(muEff,U) Gauss linear corrected; laplacian(DkEff,k) Gauss linear corrected; laplacian(DREff,R) Gauss linear corrected; laplacian(DepsilonEff,epsilon) Gauss linear corrected; laplacian((rho*(1|A(U))),p) Gauss linear corrected; laplacian(alphaEff,e) Gauss linear corrected; } interpolationSchemes { default linear; } snGradSchemes { default corrected; } fluxRequired { default no; p; }

 Thread Tools Display Modes Linear Mode

 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 OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post Wen Long Main CFD Forum 3 May 15, 2009 09:52 MĂ­Ĺˇa MilakovovĂˇ FLUENT 0 April 1, 2008 12:09 erdem NUMECA 2 August 15, 2006 15:40 kerem Main CFD Forum 0 April 24, 2006 17:04 Jonas Larsson Main CFD Forum 16 August 7, 1998 16:27

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