Supersonic Bullet Issues
 

Hi guys!

I'm trying to simulate a .50 BMG travelling at supersonic speed using Fluent, but I'm experiencing some problems. Basically I have a 2D domain (not incredibly large, I'm using the student version) and I'm simulating only a half of the bullet in steady state. These are my boundary conditions:

INLET: velocity inlet, 514.65 m/s (more or less 1.5 Mach @ 293.15 K)
OUTLET: pressure outlet (gauge pressure 0 Pa @ 293.15 K))
UP: pressure farfield (gauge pressure 0 Pa, 1.5 Mach @ 293.15 K)
DOWN: axis boundary condition (2D space axisymmetric)
BULLET: wall, no slip

As a first guess, I decided to use a k-epsilon realizable model and scalable wall functions, energy equation obviously ON, with ideal gas law and sutherland formulation for viscosity (the fluid is air). I leave some images of my mesh below.

As I was mentioning above, I'm experiencing some issues, in particular:

1 - using a Pressure-Based solver I'm able to run the case reaching a sensible result, but the continuity residuals have problems dropping down and they remain pretty high. I'm also monitoring y+ max, min and average on the wall, in addition I'm plotting the drag coefficient of the bullet; all these values remain constant at a given point of the iterative process.

2 - using a Density-Based solver, I'm not able to run the simulation at all (all the residuals diverge except for the turbulence related ones), and this makes me think of boundary conditions mistakes maybe...

With respect to point 1, I have also doubts in regards the mesh quality in one particular point (showed in two of the images below), which can represent a sort of "singularity" where the solver is not able to predict correctly the flow, since the cell should be somehow adjacent to the wall, but this does not happen due to the block construction in use (please, let me know if there is a way to do better without increase too much skewness and non-orthogonality!)

Thank you so much! :D

google for inflation layer for CFD.
you may make a small curve instead of sharp corner at front

Thanks Alexander!

Well, I tried to insert inflation layers, but it seems that Meshing does not accept them, since the mesh have already a hexa-block structure (if this is true, I think it makes sense). Obviously, I can try to adopt an unstructured mesh exploiting inflation layers at the wall (and I will try for sure), but I was wondering if the problem was solvable maintaining full hexa-block structure.

It is very common for a simulation to proceed with one solver and then diverge instantly with the other. This unimodal experience happens easily when you calibrate yourself to one type of solver and don't learn how to use the other. There is a lot more work needed than simply clicking the pressure-based or density-based check boxes, you need to have a thorough understanding of their stability characteristics as well (urfs in pressure-based solver and courant numbers in density based-solvers).


You can still use inflation layers with the structured-block hex grid, this is what results in the wicked C-type grid. However the skewness will always remain no matter your inflation layers as long as you use the same structured-block approach. You can alleviate the skewness somewhat by not using a cartesian grid and using a more curvilinear grid but still you will always have skewness in the worst possible location (at the corner) as long as you have this 3-way singularity there. If the skewness is the root of all evil for this problem then certainly invest the time to fixing it. And btw yes we do call these topological singularities because the number of blocks attached to the vertex is less than the number of edges and recognizing them will make you a good mesher.


I don't want to emphasize too much the skewed cells at the corner in case it is not the problem, but ideally with no singularity you have (in 2D) 4 blocks attached to 4 edges per vertex. The next best cases is 3-way and 5-way singularities. You have a 3-way singularity, which is not the worst. The issue is you have this singularity exactly at the corner where the Courant number is the highest because the flow needs to accelerate around the corner. Going to a unstructured-block allows you to move this singularity into the core and away from the walls and alleviates the issue by putting the bad cells in better locations.


Still, if the issue is you need to find a better initialization then don't waste all your time meshing.

Thank you LuckyTran!

You're right, I have more familiarity with the Pressure-Based solvers for sure, since I normally deal with incompressible flow simulations... I wanted to try also with the Density-Based since I know that they have been developed originally for compressible flows at high speeds, but I know also that they work in a different way with respect to the Pressure-Based ones and I need to study a little more. In particular, I have always related the Courant Number to Von Neumann stabilty analysis, dealing with transient simulations; in this case I'm performing a steady state one, and when I saw the Courant parmeters within the Fluent Control panel I remained a little bit confused :o
In any case, in regards the boundary conditions, do you think they are appropriate, or is it maybe better to work with others?


I knew that particular situation could be an issue, but I didn't know the exact definition; thanks a lot for your clarification! But at this point, how can I know if the problem is mesh or initialization related? I mean, for "initialization", are you speaking in general about all the parameters to be set (BCs, Methods, Solver etc.) or initial conditions obtained by the initialization process (by the way, I'm using the standard one, evaluating the values starting from the inlet)?

Thanks! :)

The BCs look, okay. One issue to look out for is that a velocity inlet admits both a subsonic and supersonic solution for the same velocity and that's one reason to not use it in a compressible case. But if you know this, you can still use. Just verify that you end up with the correct solution. Considering using a mass flow inlet instead. Other than this caveat, the only other BC of concern is the farfield one (because it's not a primitive BC). All the other BCs, including the velocity inlet, are primitive BCs and those are generally more robust.

Plot the solution as you solve it and the problematic areas are usually staring at you in the face. If it's a skewed cell then you see unusually high/low values of everything at that one cell and not any of its neighbors. If it's a boundary condition issue then you see sharp changes in values near the problematic boundary condition. If it's an initial condition then you see the solution struggle to evolve the right way.

The Courant number in a steady density based setup is not a typo. Even for steady cases, the solution in density-based ( and pressure-bsaed) solvers flows and there is a concept of Courant number even if it is not linked to the physical time step size. This rate of flowing is related to the radius of convergence of the solver. It should be no surprise if you simply recall that the Courant number is the speed at which information travels through the computational mesh.

Is it somehow related to the pseudo time-step used by a pseudo-transient approach or is it something more general? I know that the Courant number is related to the speed at which the information travels through the domain, but everything is evaluated with respect to the velocity at which I am advancing the computation (so, mesh characteristic length over the time-step adopted).

Substantially, a steady solver time marches (like a transient solver) to reach the steady state solution. However it does it very fast, i.e. in the case of the fluent coupled solver, with a high courant number. Also, since you do not actually care about the transient, the timestep is not global, instead it changes cell by cell so as to satisfy your imposed courant number

Thanks LoGal!

Ok, I was getting confused... When I am solving for a steady state solution, my algorithm is marching in order to reduce the imbalance within the equations; from what I know, this imbalance can be evaluated like a (pseudo) transient term. At the end of the story, using a Pressure-Based solver I get different pseudo time-steps across the cells applying fixed URFs, instead, using a Density-Based solver I end up with the same situation fixing a Courant Number. Is it correct?