CFD from the ground up

How we set up a detailed CFD base model to ensure its adaptability and applicability as a vital tool in years to come

Here I cover the process of getting our existing CAD model of last year's car into a new CFD environment. Nothing is particularly novel in this process, and the post is a bit lengthy, but there are bits and pieces throughout that will be important for and referenced by future posts. 

Why last year's and not this year's car? We need to be able to validate our simulation setup before we can trust it enough to pour in hours of design work on something new, and with this year's car currently existing only in CAD and the minds of its makers, this is our only option. 

Geometry

The first step was to simplify the car. For some parts this process is obvious; you aren't going to be needing the pedal box or the accumulator mounting bolts (for example) in an aerodynamic simulation. All these "internal" parts can simply be removed, and any mostly-enclosed hollow areas can be filled in. For other parts, the decision to include them (and if so how much to simplify them) is less obvious. Items like the wishbones, wishbone clevises, headrest, rear axle, parts housed inside the wheel structure, and even small gaps in the bodywork can be expected to affect airflow, but the question is by how much? I decided to include them all with minor simplifications to start with, as they could easily be removed from the model later if found superfluous. At the end of the simplification process, this is what I was working with (this was not the final accepted geometry after verification tests, but it does not differ too much):

Figure 1: Simplified CAD model of EV23 for CFD use with an initial guess of ideal geometric complexity

The axes and points aren't important here, they're just related to how I set up suspension motion. Notice that I've included features such as the headrest and rear suspension shocks and rockers, that clearly will have negligible effect on airflow and subsequent performance prediction of an undertray and any side-mounted aerodynamic elements that we are limiting ourselves to in the design scope for this year. This decision falls under our future planning, where we hope to carry the base CFD model forward each year, needing to change only what has also changed in reality. At some point we hope to have rear-mounted aerodynamic elements (rear wing, etc), which such geometric features will likely interact with. With a considerable amount of time planned for verification and sensitivity studies on the CFD model, inclusion of these parts was considered important, pending analysis of initial simulation results. The driver model proportions and arm positions were measured from our average-height driver sitting in the car. Given the flow through and behind the cockpit is highly disorganised, I'm not concerned about the accuracy of exact shapes or curves; I just want approximately the right flow blockage area in the right places.  

An important factor in deciding how much effort should be put into maximising CFD accuracy is that our team does not have access to a full-scale vehicle wind tunnel (or even a suitable model-scale tunnel) within reasonable means. This limits us in design complexity as it is difficult to validate complex surfaces and flow structures at a track environment, and vortex behaviour are notoriously inaccurately modelled by CFD software, particularly when multiple vortices interact. Our best option is to make the most out of what we can achieve using CFD, and incorporate this limit into our design scope. 

The first case I set up was a straight-line scenario at design ride-height. This didn't fall into our list of design cases, but having it in CFD serves as a good base model for both debugging, and importing new parts into before adjusting them to suit one of the other specific scenarios. It will also be our primary validation/correlation test case as it's the most easily reproducible driving scenario. Based on previous track data, I selected a representative speed for this case (and also the braking case as mentioned in the first post) of 18 ms-1 (65 kph). This is a good compromise between driving fast enough to provide sufficient forces (once aerodynamic components are on the car) to allow accurate measurements, while remaining slow enough to get a decent number of data points to average per run.  

The domain for the straight and braking cases is a simple rectangular prism. Dimensions were estimated erring on the larger side until domain size independence could be assessed later. Under the assumption of symmetrical flow, both these cases utilised a half-car model and symmetry plane at the centre line. For the cornering cases, this assumption was no longer valid, so the full car had to be simulated, positioned in a cylindrical domain centred at the middle of the corner (not exactly the middle though, as you'll see in the next post). 

Suspension positions for each case were taken from our custom suspension solver/simulation software, which was validated against data from several track days with last year's car.

Motion definitions

Next I had to make sure the wheels could spin properly. This is easy for tangential motion (e.g. the tyre); it's just a check-box setting. For parts with any normal motion component (e.g. the wheel centre and its "spokes"), it's a little more involved. Such parts need surrounding with a separate fluid region with its own reference frame (alternatively if running a transient simulation then the part can be made to physically rotate in the global reference frame, but for full-vehicle simulation, steady-state/RANS is the way to go for computational simplicity and so we don't have that option). With separate volumes made to specify these rotating regions (Figure 2), the region-region interfaces between the main domain and the individual moving reference frames (MRF) needed creating. The axis of rotation is shared with the tangential motion surfaces for each respective wheel. The exact position of the "spokes" would theoretically effect the result given the MRF creates only instantaneous relative motion, but it was assumed that the frequency of rotation is fast enough that any periodic effects in the wheel wake would have diffused to within any underlying uncertainty in the simulation by the time they reached critical downstream surfaces. As a retrospective note, the rear wheel MRF was later removed (replacing with the same wheel centre geometry but with no normal motion) for setup simplicity, simulation stability, and a slightly faster solve time, after the resulting flow changes were found to have negligible effect on key performance parameters. 

Figure 2: Front left wheel region encompassing all parts with a normal motion component, and assigned to a local moving reference frame (MRF)

For the cornering cases, the domain itself needed to have its own moving reference frame. There's no way to specify having "curved" velocity in the domain, so instead the domain itself needs to be rotated through the air. This was done mostly in the same way as the wheel MRF regions, except now the wheel MRFs needed to be nested within the domain MRF, which is now separate to the global reference frame (otherwise the wheels would turn on-the-spot rather than moving with the car). For both half-car and full-car cases, the walls and roof were all set to symmetry planes, which in combination with having these walls sufficiently distant to the vehicle allowed a far-field/freestream boundary approximation. 

Next up was setting surface roughness values, applied to the tyres and ground. Like everything else so far these were also a guess, consisting of a best guess of a representative sand roughness height for each. Obviously these are also subject to change later. 

Mesh

Mesh sizes were selected for all the surfaces based on a combination of research on similar applications, and a bit of intuition. These ranged from 2.5 mm for complex curvature and proximity refinements, to 5mm for wheel internals, chassis tubes, tyres, wishbones, etc, to 10 mm for bodywork, all the way up to 320 mm near the boundaries. Additional finer mesh sizes down to ~0.15 mm were then used at the tyre contact patches which had a 1 mm intersection with the ground and a 1 mm tall gap-fill around the perimeter to remove the sharp contact angle for cleaner meshing (tyres were completed changed later on - see Squishy tyres). Finally, the ground, which was otherwise set to 160mm surface sizing, was set to 5mm surface sizing in a small region directly surrounding the vehicle. 

Additional refinement volumes were added where complex (large value changes over a short distance) off-body flow was expected. These included a 20 mm sizing behind the front and rear wheels, 5 mm for the wheel internal volume, 40 mm sizing directly in front of the car, 40 mm sizing extending from the car rearwards to the domain outlet (both occupying an area a little larger than the cross section of the car) and several gradually larger volumes with gradually coarse mesh expanding out from the tyre contact patches to smoothly meet the surrounding mesh size. At this point, the half-car model had about 7 million cells, which is a big number, but quite low for a full-detail vehicle simulation. Mostly because I have no wake refinements yet other than a very crude and rather coarse rectangular volume behind each wheel. More will come later. I'll also point out that all surfaces are using a wall-modelled approach at this point with y+ >30 based on some quick research suggesting that for vehicle bodywork and non-critical geometry, wall models are perfectly sufficient.

Physics

Physics settings were generally based on default values for the k-omega SST turbulence model, with a few adjustments based on some quick research and reading of the manual:

Both values (and changes to other default values) are subject to future validation studies, particularly since a low Reynolds number, high-lift FSAE car does not represent a typical industry vehicle aerodynamics simulation scenario, from which these values were taken. Another potential alteration for racecar CFD was the addition of a no-blending option for the flux difference splitting method, however this typically only has an effect once the incompressible flow assumption begins to give way, which would be quite the feat on an FSAE car that typically runs no faster than 120kph. As such, it was assumed that this option was not required/applicable. 

Automated mesh refinement 

The final part of the setup involved creating functions for automated mesh refinement (AMR). The complexity of airflow around a vehicle (particularly with open wheels) means using simple volumetric refinements for all mesh sizing is inefficient, and cells end up being refined when they don't need to be. Some research suggested we should be looking to refine cells with high vorticity and velocity gradient. Vorticity was defined by the Lambda-2 criterion being less than 0, with additional weightings from wall distance and vorticity magnitude in Z to avoid picking up too many cells where vorticity was dominated by shear rather than well-defined stream-wise rotational structures, given this was to be a relatively fine refinement size (5mm). For all the AMR criteria, it was again a case of estimating an initial value which was later adjusted based on convergence of results monitors as part of verification testing. These later tests also led to the decision to additionally refine areas of low total pressure and high turbulent kinetic energy for decreased sensitivity of blown-up, weak, or poorly formed vortices from the suspension clevises and cockpit opening, higher visual fidelity (for easier diagnosis and hopefully easier track-test correlation), and reducing turbulence dissipation occurring too quickly. An additional buffer of 5 cells surrounding the refinement volumes were also refined to improve mesh size transition smoothness, and further reduce the positional sensitivity of high gradients (particularly for TKE and total pressure). Some of these refinement functions required a little more math and a few extra "if" statements for the cornering cases, as the velocity relative to the car differs depending on the local radius within the domain. Figure 3 shows an example refinement volume for turbulent kinetic energy, which is allowed to refine down to 5mm (give or take, depending on what the base size is set to).

Figure 3: Isosurface showing the cells flagged for refinement by one of the custom automated mesh refinement functions based on turbulent kinetic energy for a half-car simulation. The main purpose of this refinement was to reduce turbulent dissipation over Z distance, reducing what was otherwise almost immediate blending/softening of wake feature perimeters and weak vortices with freestream flow. 

Assessing results

Initial simulations were mainly used for debugging purposes, guiding small mesh adjustments to eliminate instabilities and target desired y+ values. AMR update intervals were also adjusted and set based on what we have started calling "pseudo-convergence", meaning that the simulation has converged for the current mesh refinement state, but may not be converged to the overall required tolerance. At this point, the final mesh count for a half-car case was around 14M, and 27M for a full car. I'd expect aerodynamic surfaces and the associated increase in flow complexity to add an additional 20% or so to each of these. Convergence was based on monitors for lift coefficient, drag coefficient, and aero balance being steady within 3 decimal places, and the residuals displaying no worrying characteristics. This typically took around 800 iterations, with an additional 200 iterations run for monitor and scalar displayer averaging. With 72 CPU cores, this only takes 1.5 hours for a full-detail half car case, and 3 hours for a cornering case. We have up to 128 cores available but there are significantly diminishing returns over 72 cores for our mesh count, and we're not the only ones using the shared HPC.  

The easiest way to quickly assess how the setup was performing was by using a number of plots and scenes. One of the most important checks was making sure the y+ value fell within an acceptable range on all important surfaces. Without creating hundreds of unique mesh controls it's not possible to achieve perfect y+ values for every car or ground surface, so we need a compromise. The best way to introduce and monitor this compromise is through weighting local values with the surface friction coefficient. In regions where there is no flow attachment or where velocity is low (and thus surface friction is also low), y+ can be sacrificed with negligible effect on the result. Similarly, surface area is also a good weight. It's more critical to get good y+ on a small number of very large cells than a small number of small cells, for example. I set up a series of histogram plots to show y+ weighted by the product of friction coefficient and surface area for the critical surfaces: bodywork/aerodynamic surfaces, tyres, and ground. Figure 4 shows the results after some prism layer tweaking. Some cells unavoidably remain between 5 and 30, but their weighting is small. Very high values (200+) were also avoided where possible. 

Figure 4: weighted y+ distributions for (top) bodywork and tyres, and (bottom) ground surface

To identify any surfaces with problematic y+, a scene was created to show the local value with categorical levels of An example of this is shown in Figure 5. The only area that stands out as potentially problematic in this example is the orange (very high y+) area on the sidewall of the front tyre, however tests with a locally thinner first layer height showed negligible difference in local friction, or to sidewall separation/wake shape. 

Figure 5: surface y+ values shown with categorical colour assignments: blue: y+ <5, red: 5< y+ <30, green: 30< y+ <300, orange: 300< y+

Another scene was created to check R+ (dimensionless roughness) values, and in particular the ratio of R+ to y+. Having R+ > y+ is theoretically acceptable when using the k-omega turbulence model, which does not enforce maximum roughness like k-epsilon does, however it still seems like best practise to minimise any approximations leading to inaccuracies, or sensitivities that may arise from having large areas where R+ > y+. 

The last thing I'll mention is that we found some issues with the AMR process where after refinement, the flow "solution" seemed to be having a hard time adapting to and converging on the new finer cells, with monitors and residuals floating around aimlessly for a while after each AMR update. One solution was to just wait out the convergence, but this was occasionally taking upwards of 1500 iterations. Given the issue seemed to be related to how the solution was moving through the domain, adjustments to the CFL number seemed like a logical place to look. For the half-car case, it was naturally peaking at around 150-175, and for the full-car cornering cases at 200-300. We set hard limits of 150 and 250 respectively, and this sped up convergence hugely (down to the 800 iterations mentioned earlier) while having negligible effect on the final solution. It also resulted in fewer rejected iterations, which also reduced wasted time. 

I won't go into detail on convergence/verification studies (here's a link to a post from the future), but the following parameters were all investigated using the metrics of lift coefficient, drag coefficient, aero balance, lift/drag ratio, and qualitative assessment of total pressure/wake "shapes" at critical locations to assess convergence and precision:

I'll cover post-processing visualisations and utilising CFD results in other future posts.