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

Solver for hypersonic flow in OpenFoam

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

Like Tree3Likes
  • 3 Post By sahas

Reply
 
LinkBack Thread Tools Display Modes
Old   June 18, 2012, 04:58
Default Solver for hypersonic flow in OpenFoam
  #1
Member
 
Charles
Join Date: May 2011
Location: France
Posts: 78
Rep Power: 6
Santos-Dumont is on a distinguished road
Hi,

I'm trying to simulate the wake behind an object at hypersonic speeds (around M=15) with OpenFoam.
The solver SonicFoam seems inappropriate for hypersonic speeds because the shock is too close to the object. So I used rhoCentralFoam which gives nice results.

The problem is that I have a very huge domain as I'm interested in the far wake behind an object and rhoCentralFoam is a transient solver. My first case (half of a sphere with a domain lenght of 2000 times the sphere diameter) needed 10 days of calculations on 7 CPU cores. Even if the results are nice the calculation time is too important.

Which solver would you recommend for such a case ? Is there a steady-state solver which is fine for hypersonic flows ?

Thanks,
Charles
Santos-Dumont is offline   Reply With Quote

Old   June 19, 2012, 02:56
Default
  #2
Member
 
Join Date: May 2009
Posts: 31
Rep Power: 8
KrisT is on a distinguished road
I don't think there are too many options here:

In principle you could save time by using adaptive mesh (maybe you already do?), but when I last tested (on rhoCentralFoam), the capability in OpenFOAM it couldn't handle multi-core simulations in an efficient manner, i.e., there was no redistribution of cells after refinement. Which resulted in that - when using several cores - one didn't gain much (on single core it was working well though).

I would be very interested to hear if there has been improvements in OpenFOAM in this respect (did my tests on 1.5 or 1.6).

Apart from this, try to find more cores ;-)

/jon
KrisT is offline   Reply With Quote

Old   June 19, 2012, 03:26
Default
  #3
Member
 
Charles
Join Date: May 2011
Location: France
Posts: 78
Rep Power: 6
Santos-Dumont is on a distinguished road
Hi Kris,

I think I could try coarsening the mesh or adapt it but I might reduce the precision of the results as I am still trying to get a mesh convergence on that case (using another CFD software).

I had no particular trouble running it in parallel. I didn't make a measurement on the actual speed gained by the parallel run but I tried to run it on a single core and it was way, way much slower.

Charles
Santos-Dumont is offline   Reply With Quote

Old   June 19, 2012, 04:34
Default
  #4
Senior Member
 
Olivier
Join Date: Jun 2009
Location: France, grenoble
Posts: 235
Rep Power: 9
olivierG is on a distinguished road
hello,

you may try to use coEuler / SLTS time scheme, then your time step will be higer (at least ~ x5 time faster ).

regards,
olivier
olivierG is offline   Reply With Quote

Old   June 20, 2012, 02:54
Default
  #5
Member
 
Join Date: May 2009
Posts: 31
Rep Power: 8
KrisT is on a distinguished road
Hi! Sorry, but I meant you should try the _adaptive_ mesh refinement, so that mesh is updated during the calculation (according to some measure - involving gradient of density or so). That is, mesh is refined around shocks, and more coarse where not much happens in the field. But as I said, not sure if this is efficient when running in parallel.
KrisT is offline   Reply With Quote

Old   July 23, 2012, 07:21
Default
  #6
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
If you are interesting in local time stepping (use "pseudo-time") there is some information in the forum about its implementation in rhoCentralFoam. I have done it for rhoCentraFoam2.0, and I can post here the code if demanded :-)
sahas is offline   Reply With Quote

Old   July 24, 2012, 02:52
Default
  #7
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10
vkrastev is on a distinguished road
Hi Alexander,
I'm interested in rhoCentralFoam+pseudo-time stepping for accelerating convergence to steady state. Can you please post your code?

Thanks

V.
vkrastev is offline   Reply With Quote

Old   July 25, 2012, 07:16
Default
  #8
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
There are three files in attachment (modified rhoCentralFoam.C, createFields.H, new setrDeltaT.H).
Below I post diff's

rhoCentralFoam.C
Code:
--- ../../OpenFOAM-2.0.0/applications/solvers/compressible/rhoCentralFoam/rhoCentralFoam.C<---->2011-06-16 16:48:41.000000000 +0400
+++ rhoCentralFoam.C<-->2012-07-25 18:46:33.000000000 +0400
@@ -35,6 +35,7 @@
 #include "turbulenceModel.H"
 #include "zeroGradientFvPatchFields.H"
 #include "fixedRhoFvPatchScalarField.H"
+#include "fvcSmooth.H"
.
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
.
@@ -149,8 +150,9 @@
         // estimated by the central scheme
         amaxSf = max(mag(aphiv_pos), mag(aphiv_neg));
.
-        #include "compressibleCourantNo.H"
         #include "readTimeControls.H"
+        #include "setrDeltaT.H"
+        #include "compressibleCourantNo.H"
         #include "setDeltaT.H"
.
         runTime++;
createFields.H
Code:
--- ../../OpenFOAM-2.0.0/applications/solvers/compressible/rhoCentralFoam/createFields.H<------>2011-06-16 16:48:41.000000000 +0400
+++ createFields.H<---->2012-07-25 18:48:04.000000000 +0400
@@ -111,3 +111,19 @@
         thermo
     )
 );
+
+volScalarField rDeltaT
+(
+    IOobject
+    (
+        "rDeltaT",
+        runTime.timeName(),
+        mesh,
+        IOobject::NO_READ,
+        IOobject::AUTO_WRITE
+    ),
+    mesh,
+    1/dimensionedScalar("maxDeltaT", dimTime, GREAT),
+    zeroGradientFvPatchScalarField::typeName
+);
+
setrDeltaT.H
Code:
--- ../../OpenFOAM-2.0.0/applications/solvers/compressible/rhoCentralFoam/setrDeltaT.H<>1970-01-01 03:00:00.000000000 +0300
+++ setrDeltaT.H<------>2012-07-25 18:50:06.000000000 +0400
@@ -0,0 +1,28 @@
+{
+
+    // Set the reciprocal time-step from the local Courant number
+    rDeltaT.dimensionedInternalField() = max
+    (
+      1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT),
+      fvc::surfaceSum(amaxSf)().dimensionedInternalField()
+              /(maxCo*mesh.V())
+    );
+
+    Info<< "    Flow        = "
+            << gMin(1/rDeltaT.internalField()) << ", "
+            << gMax(1/rDeltaT.internalField()) << endl;
+
+    // Limit the largest time scale
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    rDeltaT.max(1/maxDeltaT);
+
+
+    // Spatially smooth the time scale field
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    scalar rDeltaTSmoothingCoeff = 0.1;
+    fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
+
+    Info<< "    Overall     = " << min(1/rDeltaT).value()
+        << ", " << max(1/rDeltaT).value() << nl << endl;
+
+}
To activate Local Time Stepping (LTS) you should specify the next option in system/fvSchemes:
Code:
ddtSchemes
{
    default         localEuler rDeltaT;
}
Apropos if you need the "waveTransmissive"boundary condition with LTS, you should correct the file OpenFOAM-2.0.0/src/finiteVolume/fields/fvPatchFields/derived/advective/advectiveFvPatchField.C in the next way
Code:
--- ../OpenFOAM-2.0.0/src/finiteVolume/fields/fvPatchFields/derived/advective/advectiveFvPatchField.C    2011-06-16 16:48:41.000000000 +0400
+++ /opt/openfoam200/src/finiteVolume/fields/fvPatchFields/derived/advective/advectiveFvPatchField.C    2012-04-08 21:32:23.000000000 +0400
@@ -30,6 +30,7 @@
 #include "EulerDdtScheme.H"
 #include "CrankNicholsonDdtScheme.H"
 #include "backwardDdtScheme.H"
+#include "localEulerDdtScheme.H"
 
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
@@ -226,6 +227,7 @@
         (
             ddtScheme == fv::EulerDdtScheme<scalar>::typeName
          || ddtScheme == fv::CrankNicholsonDdtScheme<scalar>::typeName
+         || ddtScheme == fv::localEulerDdtScheme<scalar>::typeName
         )
         {
             this->refValue() =
Also in the version 2.0 there is a bug with Prandl number, which is alway 1 (see bug http://www.openfoam.org/mantisbt/view.php?id=475)
Attached Files
File Type: c rhoCentralFoam.C (7.6 KB, 37 views)
File Type: h createFields.H (2.1 KB, 28 views)
File Type: h setrDeltaT.H (825 Bytes, 31 views)
sahas is offline   Reply With Quote

Old   July 25, 2012, 08:26
Default
  #9
Senior Member
 
Vesselin Krastev
Join Date: Jan 2010
Location: University of Tor Vergata, Rome
Posts: 361
Rep Power: 10
vkrastev is on a distinguished road
Thanks a lot, I'll give a try to the modified code as soon as possible.

Regards

V.
vkrastev is offline   Reply With Quote

Old   August 27, 2012, 09:05
Default
  #10
ndr
New Member
 
Join Date: Aug 2009
Location: Stuttgart, Germany
Posts: 20
Rep Power: 7
ndr is on a distinguished road
Hi Alexander,

I originally had a slightly different solution for the local time-stepping implementation in rhoCentralFoam for OF 2.1.x, only my rDeltaT was calculated a bit different.

However, when I tried to validate the solver for a simple flat plate case with a laminar Mach 2 flow, the results were not completely satisfying. So I gave your solver a try and interestingly nearly the same behaviour occurred: Especially the pressure field is fluctuating wildly downstream of the shock wave which originates at the plate's leading edge. I attached the contour plots of the pressure field both for the case using LTS and for exactly the same setup using the normal timestepping. Without LTS only one shock wave is captured as it should be. For the boundary conditions: The flow is coming from the left with uniform velocity, pressure and temperature. I used wave transmissive for the right outlet and the top, the plate at the bottom starts at 50 cm (before that a symmetry plane is used). For the reconstruction schemes I used the Gamma scheme.

As I had used the original setup for validation of the boundary layer I also had a look into that. And as you can see in the velocity and temperature profiles across the BL at a length of 1 m (red is without LTS, blue is with LTS) those profiles are not really identical... The red profiles are validated against literature and an analytical solution, so they are definitely correct.
I already tried to change the time-step smoothing, but that had no effect at all.

Do you have an idea where this fluctuation may come from? Or did you maybe observe similar problems when testing your solver?

Regards

Nils
Attached Images
File Type: jpg LTS-off.jpg (13.0 KB, 98 views)
File Type: jpg LTS-on.jpg (13.6 KB, 85 views)
File Type: jpg LTS-UProfile.jpg (24.5 KB, 63 views)
File Type: jpg LTS-TProfile.jpg (25.6 KB, 46 views)

Last edited by ndr; August 28, 2012 at 02:20.
ndr is offline   Reply With Quote

Old   August 28, 2012, 05:25
Default
  #11
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
Hi Nils,

what Courant number "maxCo" do you use? I found that sometimes there are fluctuations in the fields when Courant number is greater than approximately 0.1.
sahas is offline   Reply With Quote

Old   August 28, 2012, 08:24
Default
  #12
ndr
New Member
 
Join Date: Aug 2009
Location: Stuttgart, Germany
Posts: 20
Rep Power: 7
ndr is on a distinguished road
Hi Alexander,

ok, I didn't realize I had to take a Courant number that low. I ran my calculations as usual with a maxCo of 0.4, but I also tried 0.2 without any improvement in the results.
But as you suggested I will also try 0.1 (and 0.05) and will let you know if that helps with the fluctuations.
ndr is offline   Reply With Quote

Old   August 28, 2012, 09:21
Default
  #13
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
Also, the result may strongly depend on interpolation scheme used. For example, vanLeer sometimes gives oscillations, whereas upwind does not. As well you may try SFCD scheme.
P.S. Of course, interpolation scheme is weakly related to ddt scheme, but who knows exactly?..
sahas is offline   Reply With Quote

Old   August 28, 2012, 10:43
Default
  #14
ndr
New Member
 
Join Date: Aug 2009
Location: Stuttgart, Germany
Posts: 20
Rep Power: 7
ndr is on a distinguished road
Yes, I already noticed that vanLeer seems to be unstable for certain problems, especially when dealing with hypersonic boundary layers.
But until now the Gamma scheme always worked fine for me, so I think I will try to keep as many parameters constant as possible (when comparing localEuler and Euler) so that differences in the results can easily be traced back to the ddt scheme.
ndr is offline   Reply With Quote

Old   August 29, 2012, 08:22
Default
  #15
ndr
New Member
 
Join Date: Aug 2009
Location: Stuttgart, Germany
Posts: 20
Rep Power: 7
ndr is on a distinguished road
Hi Alexander,

thanks for the hint regarding the maximum Courant number, that did the trick. If I use maxCo 0.1 or lower it works perfectly well.
It's a bit of a drawback that the number of iterations now increases considerably compared to the original Euler case, which was calculated at maxCo = 0.4 and without any problems. But the LTS still seems to converge faster to a steady state by a factor of about 2.
ndr is offline   Reply With Quote

Old   August 29, 2012, 11:06
Default
  #16
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
Glad to hear good news, Nils

As I understand, the "drawback" is a kind of "phantom". Look: when you use Euler scheme you limit the time step corresponding to Courant number 0.4 in some cell (may be or greatest, or smallest, or anything else). But in another cell Courant number corresponding to chosen time step may be much less than 0.4 (for example, 0.04) - and this cell may be source of oscillations. Whenever you use localEuler, in all cells (almost) you have Courant number 0.4, also in "problem" cell, and the critical Courant number for it cannot be greater than 0.1. This is just my theory, I did not have "hard look" on it
sahas is offline   Reply With Quote

Old   August 29, 2012, 14:13
Default
  #17
ndr
New Member
 
Join Date: Aug 2009
Location: Stuttgart, Germany
Posts: 20
Rep Power: 7
ndr is on a distinguished road
Yes, I was aware of that. However, due to my grid structure the size of the cells in flow direction (which should be mainly responsible for Courant number calculation as there is hardly any crossflow) was approximately the same all over the grid, only the wall normal distance was refined further.

So when I ran the Euler calculation I indeed did not get a Courant of 0.4 all over the grid, but the discrepancy between the cells was not that big. If I remember correctly the minimum Courant number was between 0.05 and 0.1, but that occured in a region far away from the wall and any shocks, so those cells should not be critical for any oscillations.

But it seems that somewhere around Co = 0.1 there is the upper limit for LTS, as I tested my old version of the solver again and it showed the same behaviour. The speed-up is still considerably high though and for quasi steady-state calculations it seems to be the only reasonable choice for rhoCentralFoam, all other simulations may take ages to compute depending on your geometry and flow type....
ndr is offline   Reply With Quote

Old   August 30, 2012, 03:20
Default
  #18
Member
 
Alexander
Join Date: Mar 2009
Posts: 39
Rep Power: 8
sahas is on a distinguished road
By the way, inside boundary layer the diffusive term has significant value so I think one should consider also the diffusion Courant number (mu*dt/dx^2) - may be, it is responsible for arising oscillations?
sahas is offline   Reply With Quote

Old   August 30, 2012, 09:31
Default
  #19
New Member
 
S. Mitra
Join Date: Oct 2011
Location: US
Posts: 1
Rep Power: 0
SMitra is on a distinguished road
Hi,

You may try to use normalized properties for this compressible problem. That way, you may be able to increase the time step that you need to solve the problem. This will give you quicker qualitative results for your problem.
For correct quantitative solution, you will have to use the proper properties.

Hope this helps.

PS: For normalized prop: look at "thermophysicalproperties" in the ForwardStep problem in sonicFoam tutorial.
SMitra is offline   Reply With Quote

Old   October 21, 2012, 01:42
Default
  #20
Super Moderator
 
praveen's Avatar
 
Praveen. C
Join Date: Mar 2009
Location: Bangalore
Posts: 244
Blog Entries: 6
Rep Power: 9
praveen is on a distinguished road
Quote:
Originally Posted by sahas View Post
By the way, inside boundary layer the diffusive term has significant value so I think one should consider also the diffusion Courant number (mu*dt/dx^2) - may be, it is responsible for arising oscillations?
The diffusive terms in momentum and energy equation are treated implicitly. So diffusive time step should not be a problem. Is the mesh made of quadrilaterals or triangles ?
praveen is offline   Reply With Quote

Reply

Tags
hypersonic wake

Thread Tools
Display Modes

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


Similar Threads
Thread Thread Starter Forum Replies Last Post
Solver for an incompressible, turbulent flow with heat transfer tH3f0rC3 OpenFOAM Running, Solving & CFD 6 March 24, 2011 06:41
nano fluif flow in openfoam anijdon OpenFOAM 0 February 18, 2011 17:22
What solver to use for inviscid flow simulation over missile mecbe2002 OpenFOAM 0 April 27, 2010 11:10
What solver to use for simulating flow behavior in inlet manifold mecbe2002 OpenFOAM 4 April 9, 2010 08:01
OpenFOAM Training in Europe and USA hjasak OpenFOAM 0 August 8, 2008 05:33


All times are GMT -4. The time now is 13:45.