Publications MScThesis Carreira 2022
Publications MScThesis Carreira 2022
Publications MScThesis Carreira 2022
Aerospace Engineering
Examination Committee
Chairperson: Prof. Filipe Szolnoky Ramos Pinto Cunha
Supervisor: Prof. André Calado Marta
Member of the Committee: Prof. José Manuel Da Silva Chaves Ribeiro Pereira
June 2022
ii
Dedicated to my parents
iii
iv
Acknowledgments
The development of this dissertation would not be possible without the help and support of many peo-
ple. I would like to express my gratitude and appreciation for both my supervisors Professor André Marta
and Professor Luı́s Eça whose guidance, analysis, insights and encouragement have been invaluable
throughout this study.
Many thanks to the Formula Student Técnico team for making this work possible. Much effort has
been made by team members to organize and execute a testing day for aerodynamic validation and I am
grateful to them all. Particular thanks go to André Correia, Pedro Oliveira, Miguel Valverde and Diogo
Jerónimo for integration with the current team and test planning, Bernardo Taveira for helping with data
handling from the test and João Perdigão for assisting in the preparation and development of a plan
for aerodynamic testing. I would also like to thank my friends Jaime Pacheco and José Luciano for our
time in FST team that culminated in FST10e aerodynamics development, which enabled me to continue
further exploring it during this work.
I also wish to thank all my closest friends for their encouragement and for all throughout this academic
journey.
Finally, I would like to thank all my family, especially my parents and brother for their continued
support and trust throughout all these university years.
v
vi
Resumo
As forças aerodinâmicas de um carro de competição variam consoante a sua atitude, que por sua
vez muda ao longo da pista devido às acelerações a que é sujeito. O principal objectivo deste trabalho
é perceber como estas mudanças de atitude influenciam a aerodinâmica do FST10e, o último protótipo
construı́do pela equipa FST Lisboa.
A atitude do carro foi decomposta em 5 parâmetros de interesse: distância ao solo dos eixos di-
anteiro e traseiro e ângulos de rolamento, direção e guinada, onde o comportamento do veı́culo é
estimado por simulações CFD. O foco principal de CFD foi simular em condições de curva, em que
foi estimado o respetivo erro numérico ao realizar um estudo de convergência de malha. A influência
de cada um destes parâmetros foi analisada individualmente revelando alguns resultados inesperados,
sendo que o ângulo de guinada apresentou com a maior sensibilidade com uma mudança em 16% da
força descendente ao longo do seu intervalo. O ângulo de rolamento revelou ter um comportamento
monotónico e algo peculiar resultando numa alteração de 9% na força descendente e uma deslocação
do centro de pressão de 12%. A distância ao solo apresentou um mudança significativa de 10% em
CD .A e 20% em −CL .A entre valores máximos em condições simétricas. Um mapa aerodinâmico foi
então criado tendo por base mais de 100 pontos de dados contribuindo para uma melhor compreensão
sobre a janela de funcionamento. Os resultados podem ser obtidos rapidamente usando um modelo
aproximado, que supera o processo moroso de análise CFD. Por fim, uma componente de validação
aerodinâmica também foi abordada para estimar os coeficientes de força descendente e de arrasto com
testes em pista. O maior desafio encontrado foi estimar a velocidade do ar durante os testes, que se
revelou uma grande fonte de incerteza para correlacionar os resultados em pista com CFD.
vii
viii
Abstract
Aerodynamic forces on a race car depend on its attitude, which changes along the track given the
accelerations they are subjected to. The main objective of this work is to understand how these attitude
changes affect the FST10e aerodynamics, the latest prototype made by FST Lisboa team.
The car’s attitude was broken down into five parameters of interest: front and rear ride heights, roll,
steering and yaw angles, with the vehicle behaviour being estimated by CFD simulations. The primary
focus of CFD was cornering condition simulation, with the corresponding numerical error estimated by
performing a mesh convergence. The influence of each of these parameters was evaluated individually,
revealing some unexpected results, with the most impactful sensitivity being the yaw angle with a 16%
change in downforce across the interval. Roll angle revealed to have a monotonic and rather indepen-
dent behaviour resulting in 9% change in downforce and a shift of 12% in pressure centre. Ride height
presented a significant change with a 10% change in CD .A and 20% in −CL .A between peak values
in symmetric conditions. An aerodynamic map was then created based on over 100 data points, which
provided a thorough understanding across the whole working envelope. The results can be promptly
obtained by using a surrogate model, which overcomes the slow CFD analyses. Lastly, an aerodynamic
validation component was also addressed to estimate downforce and drag coefficients with on-track
testing. The main challenge encountered was to estimate air speed during validation testing, which
revealed to be a great source of uncertainty to correlate CFD with track results.
Keywords: CFD, Formula Student, Aerodynamic mapping, On-track validation, Ride height,
Vehicle attitude
ix
x
Contents
Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . v
Resumo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii
Abstract . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix
List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xv
List of Figures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xvii
Nomenclature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xxi
Glossary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xxiii
1 Introduction 1
1.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2 Topic Overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2.1 Race Car Aerodynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.2.2 Formula Student . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.2.3 Aerodynamic Mapping . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4
1.3 Objectives and Deliverables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.4 Thesis Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2 Mathematical Modelling 9
2.1 Geometric Modelling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
2.1.1 Coordinate System . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
2.1.2 Geometry . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.1.3 Parameters of Interest . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2 Cornering Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.2.1 Pre-Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2.2.2 Solver . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
2.3 Adaptive Mesh Refinement . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
2.4 Numerical Error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.4.1 Round-Off Error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.4.2 Statistical Sampling Error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
2.4.3 Iterative Error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.4.4 Discretization Error . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
xi
3 Simulation Setup 33
3.1 Straight Line Simulation Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
3.1.1 Lateral Symmetric Conditions Using Half Vehicle . . . . . . . . . . . . . . . . . . . 33
3.1.2 General Conditions Using Complete Vehicle . . . . . . . . . . . . . . . . . . . . . . 34
3.2 Process Automation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
3.2.1 Post-Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
3.2.2 Macro . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
4 Parametric Analysis 43
4.1 Ride Height . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
4.2 Roll Angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
4.3 Steering Angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
4.4 Yaw Angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
4.5 Results Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
5 Surrogate Model 51
5.1 Design of Experiments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
5.2 Construction of Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
5.3 Validation of Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
5.4 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 55
6 On-track Validation 59
6.1 Testing Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
6.1.1 Constant Speed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
6.1.2 Coast Down . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
6.1.3 Instrumentation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
6.1.4 Testing Procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
6.2 Data Processing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
6.2.1 Loaded Tyre Radius . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
6.2.2 Dynamic Ride Height . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
6.2.3 Static Ride Height . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67
6.2.4 Atmospheric Conditions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
6.3 Experimental Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
6.3.1 Constant Speed . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
6.3.2 Coast Down . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
7 Conclusions 75
7.1 Achievements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75
7.2 Future Work . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
Bibliography 77
xii
A Surrogate Model Data 81
A.1 Data Points . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81
A.2 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
A.3 Parametric Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
A.4 Yaw Angle Evolution of −CL .A . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
B Validation Results 89
B.1 Constant Speed Test Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 89
B.2 Additional Information . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
B.3 Coast Down Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
B.4 Example Constant Speed and Coast Down Runs Data . . . . . . . . . . . . . . . . . . . . 90
B.5 Suspension Kinematics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97
xiii
xiv
List of Tables
xv
xvi
List of Figures
1.1 Effect of aerodynamic downforce on the polar diagram for a tire maximum performance [2] 2
1.2 FST10e components . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.3 FST10e in FSG 2021 during acceleration event . . . . . . . . . . . . . . . . . . . . . . . . 4
1.4 Formula Student track layout . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.5 Ride height envelope for a formula car [6] . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.6 Aerodynamic maps for a formula racecar [6] . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.7 Aerodynamic mapping and analysis by Total Sim [8] . . . . . . . . . . . . . . . . . . . . . 7
xvii
2.24 Downforce standard deviation and error for the last 500 iterations . . . . . . . . . . . . . . 29
2.25 Comparison of the effect of averaging in the total pressure coefficient . . . . . . . . . . . 29
2.26 Numerical uncertainty and error estimate for low and high-speed corner attitudes . . . . . 32
6.1 Military and technical training center of the air force, Google Maps . . . . . . . . . . . . . 59
6.2 Least-squares approximation of experimental data (107 V-D points): d(v) = 0.135 +
0.00104v + 0.00025v 2 (CD = 0.342) [43] . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
6.3 Xsens MTi 600-series [44] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
6.4 ETAS ES910 [46] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
6.5 Linear potentiometer . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
6.6 Tyre radius as function of vertical load [48] . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
6.7 Tyre data handling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
6.8 Workflow for potentiometer data-processing . . . . . . . . . . . . . . . . . . . . . . . . . . 66
xviii
6.9 Front ride height configurations with different shim size . . . . . . . . . . . . . . . . . . . . 67
6.10 Suspension diagram by FST Lisboa vehicle dynamics department . . . . . . . . . . . . . 68
6.11 Wool tufts during testing day . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
6.12 Overview of a testing run . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
6.13 Constant speed test results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72
6.14 Overview of coast down test data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73
6.15 Acceleration vs velocity from coast down test . . . . . . . . . . . . . . . . . . . . . . . . . 73
6.16 Coast down results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
xix
xx
Nomenclature
Greek symbols
β Angle of side-slip.
δ Angle of steering.
Ω Angular velocity.
φ Angle of roll.
ψ Angle of yaw.
ρ Density.
θ Angle of pitch.
Roman symbols
A Reference area.
b Reference length.
CD Coefficient of drag.
CL Coefficient of lift.
CM Coefficient of moment.
CY Coefficient of sideforce.
p Pressure.
q Dynamic pressure.
V Velocity magnitude.
xxi
Subscripts
∞ Free-stream condition.
I Inertial component.
i, j, k Computational indexes.
R Rotating component.
x, y, z Cartesian components.
Superscripts
b Approximate model.
T Transpose.
xxii
Glossary
xxiii
RANS Reynolds-Averaged Navier–Stokes Equations
RCH Roll Centre Height
RC Roll Centre
RHS Right Hand Side
RRH Rear Ride Height.
SSL Standard Sea Level
SST Shear Stress Transport
SS Split Sample
SoC State of Charge
TSTP Test Specification and Test Procedure
xxiv
Chapter 1
Introduction
1.1 Motivation
The aim of race cars is ultimately to win races and to do so, there are many areas where it is possible
to reduce lap times and beat the competitors. One area of expertise that has been on a rise since 1970s
is aerodynamics and a lot of performance can be extracted from it. Aerodynamics plays an important
role in the vehicle’s stability and handling, helping the driver to achieve consistent and quick lap times.
In contrast to aircraft, negative lift (or downforce) is beneficial for the car’s performance on track. For
this reason, the same principles of aircraft aerodynamics are used but instead they are turned upside
down in order to generate downforce.
As the car is in constant motion due to the accelerations it is subjected to (laterally when cornering
and longitudinally when accelerating/braking) vehicle’s attitude varies a lot across the track. These
changes in the vehicle state directly affect the car aerodynamic performance and the motivation to this
work was to understand how the aerodynamic forces and moments acting on a race car vary over a
range of attitudes.
The comprehension of this behaviour helps to understand in which conditions aerodynamics package
can be refined but is also advantageous for the vehicle dynamics point of view. It provides information to
help identify zones of interest to ride the car in and adjust car setup to ride on those conditions, reaching
the full potential of the prototype. Lastly, it is also possible to draw conclusions about the stability and
”feeling” of the car that can help identify reasons for understeer or oversteer 1 .
The only components of the vehicle that are in contact with the road are tires, more specifically the
contact patch area. All forces experienced by the vehicle are in the end supported by the tires. The
1 Terms used to describe the sensitivity to steering the car. It undergoes understeer when it does not turn enough as com-
manded by the driver and oversteer when the car turns more than desired
1
aerodynamic loads change significantly the forces in them, and they play a very important role in race
car performance.
There are three major effects of aerodynamic forces [1] :
• Downforce and its distribution: Aerodynamic vertical force enables an expansion in the tire working
envelop, as seen in Fig. 1.1, increasing the size of the polar or G-G diagram 2 . Having this area
increase, it affects maximum cornering, braking and traction accelerations. Its distribution between
front and rear wheels affects stability due to the grip difference between these wheels.
• Drag: Resistance force has an effect in maximum speed achieved by the vehicle and also plays a
part in longitudinal accelerations.
• Side force: Side component of aerodynamic forces point to the corner centre affecting centrifugal
force experienced by the vehicle, that is once again related to stability. Also associated with the car
handling and control are all moments created such as yaw, roll and pitch (another way to represent
downforce distribution).
Figure 1.1: Effect of aerodynamic downforce on the polar diagram for a tire maximum performance [2]
Aerodynamic forces are generated from distributed loads in different components across the car. In
race cars, particularly in formulas, some typical aerodynamic devices are front and rear wing, sidepods
and a floor/diffuser. Figure 1.2 describes the different components of the aerodynamic package of
FST10e prototype.
Formula Student (FS) [3] is an engineering competition for university students around the world, that
consists of conceiving, designing, fabricating, developing and competing with a small formula style race
car. The competition is divided into 3 classes:
2 G-G diagram is a way to represent the performance envelop of a car where it shows the zone in which the car remains stable
2
Figure 1.2: FST10e components
Formula Student Técnico Lisboa (FST Lisboa) recently competed with an EV and a DV-EV vehicles
following Formula Student Germany (FSG) rules. The EV competition is split into Static and Dynamic
Events [4] that include:
• Static Events
– Business Plan Presentation: develop a business model around the prototype race car to
generate profit;
– Cost and Manufacturing: understanding of manufacturing processes and cost associated with
a prototype race car but also the difference to a mass production;
– Engineering Design: evaluate the engineering process involved in the vehicle design.
• Dynamic Events
This Thesis work fits into 1 static event that is Engineering Design. In this event, a panel of judges
evaluate the engineering process and effort made by the students to design the whole vehicle. Regard-
ing aerodynamics, an aerodynamic mapping provides information about the performance in different
attitudes. It is important to design an aerodynamic package that works as intended in all working en-
velope, and integrating aerodynamic mappings into the design process helps in achieving this goal (for
example, to prevent a diffuser stalling when the ground clearance reduces and cause an unstable be-
haviour). With respect to vehicle dynamics design, aerodynamic mappings can play an important role
into lap time simulators, as they can be an input for each time step calculation. By updating the aerody-
namic behaviour along the track, it enables a better estimate of the car response and with these results,
it is possible to evaluate lap time change for multiple designs and even try different car tunings to ride
3
in the zones of interest. This is yet to be done in FST simulators, but developments are currently being
made to implement this tool in the future.
Regarding dynamic events, static ride heights and wings angles/position can be adjusted based
on the information provided by aerodynamic mappings as seen in Fig. 1.3 with the rear wing angles
adjusted for a low drag configuration for acceleration event.
A track layout of a FS competition is displayed in Fig. 1.4 together with velocity and power usage
prediction from a point mass simulator by FST’s vehicle dynamics department. It is clear that there are
only a few sections that the car runs in a straight line or is in constant turning and acceleration. These
track characteristics support the idea of studying the vehicle in all attitudes due to being in constant
motion.
Aerodynamic mapping is a visual representation of aerodynamic performance for a given set of in-
dependent variables (at least two). In motorsport, this tool is used to understand how the car behaves
across its working envelope and usually they include front and rear ride height (which represent the dis-
tance between the ground to the front and rear axis respectively). An example of a ride height envelope
is shown in Fig. 1.5 . At the bottom left corner, we have brake dip that happens when the car pitches for-
ward due to longitudinally deceleration resulting in a decrease of front and increase of rear ride heights.
4
(a) Velocity across a lap (b) Power usage over a lap
After this transitory effect, the remaining braking will progressively reduce speed and consequently aero-
dynamic forces, translating into both front and rear ride heights to increase (suspension springs extend
with lower applied forces). When the driver releases the brake pedal and starts throttling at corner entry
the opposite occurs, vehicle pitches backwards (front goes up and back goes down in top right corner)
and when it gets into to a straight line again, the speed increase compresses the suspension.
The aerodynamic performance is usually evaluated by at least three dependent variables: lift coeffi-
cient, drag coefficient and aerodynamic balance (aero balance). Lift coefficient, defined by
L
CL = 1 2 , (1.1)
2 ρAV∞
is a dimensionless number that indicates the force perpendicular to the incoming airflow. It quantifies
the downforce level of the car that change the vertical load tyres are subjected to. Drag coefficient,
D
CD = 1 2 , (1.2)
2 ρAV∞
represents the force in the direction of incoming airflow and it is directly related to the power used to over-
come air resistance and longitudinal behaviour. The variables necessary to calculate these coefficients
are: lift force L, drag force D, air density ρ, reference area A and freestream air velocity V∞ .
Lastly, aerodynamic balance is a common way to represent pressure centre location in race cars and
5
it is simply a percentage of how much front axle downforce there is over total downforce expressed by
This balance is associated with stability once the magnitude of front wheel downforce can translate into
understeer or oversteer.
An example of 2D aerodynamic maps is shown in Fig. 1.6 with CL and CD represented on the bottom
two images and aero balance is pictured in a different way by using front and rear downforce separately.
A more in depth analysis of ground effect aerodynamics of race cars can be found in [7] that display
the change in downforce in an open wheel race car with absolute values comparing to wind tunnel test
experiments.
Although using two variables can provide plenty of information, for a more thorough analysis the
number of independent variables can be increased, which is the case of this work having a total of five
that will be described into detail during Section 2.1.3. In Fig. 1.7(a) there is an example of how can
an aerodynamic mapping look like when it includes more than two independent variables and, in Fig.
1.7(b), it is also possible to understand how can downforce and drag vary across an entire map.
Usually these complex mappings can be approximated by a surrogate model
The main objective of this work is to create an aerodynamic map of the FST10e prototype, the latest
vehicle created by FST that attended competitions in 2021. It is intended to include a maximum of five
parameters - front and rear ride heights, roll, steer and yaw angles - and use full car CFD simulations.
6
(a) Contour plots showing changes in performance (b) A plot of drag vs. downforce of an entire map
Another major objective is to do on-track testing and include a validation component. This step
would bridge the gap between CFD results and reality, allowing to obtain a comparison between the
two. The goal with these track tests is to primarily focus on trends between different setups and not on
absolute values correlation due to the uncertainties that cannot be accounted for and also some limited
resources to achieve it. Nonetheless, trend evaluation would be a major step towards validation focus
which is currently one of the key objectives for FST Lisboa.
In the course of this Thesis, it is expected to present the following deliverables:
• Parametric CAD: Making a parametric geometry with respect to the variables in study to speed
up the process of making multiple iterations;
• CFD simulation: Setup three types of simulations that include cornering, half car symmetrical and
full car in straight line giving more emphasis to corner type;
• Mesh convergence: Estimate numerical error from a mesh convergence study;
• CFD process automation: Develop a macro to fully automate CFD simulations to reduce repeti-
tive work hours and minimize possible human input errors;
• Aerodynamic map: Present an aerodynamic mapping with a reduced number of simulations, as
they are computationally expensive;
• Correlation CFD-track results: Perform on-track testing and compare it to simulation results.
During this Chapter 1, a brief introduction to the work was made. Initially, it was explained the mo-
tivation that led to this work and then some background was shown regarding race car aerodynamics,
Formula student and aerodynamic maps. At the end of the present Chapter, the objectives and deliver-
ables for this work were set and the outline of this Thesis.
7
Next, there is Chapter 2 where the mathematical model is described. It starts with the car geometry
followed by how corner simulation were modelled and respective numerical error estimation.
In Chapter 3 the remaining type of simulations setup are presented together with the CFD automation
workflow.
In Chapter 4 it is assessed the individual impact of each variable of interest with parametric studies.
Chapter 5 focus on the Surrogate Model explaining the approach used to build an approximate model
based on CFD results in order to reduce computational cost.
Chapter 6 addresses the on-track validation tests, detailing the testing procedure, data processing
gathered from the tests and the results obtained.
Finally, Chapter 7 is the conclusion that states the achievements in this work and also future work
that can be done to improve and expand knowledge in this subject.
8
Chapter 2
Mathematical Modelling
In this Chapter, the models used to obtain a numerical solution for this work will be described from
geometric modelling to the corner simulation settings and at the end a numerical error estimation will be
addressed.
To study the effects of aerodynamics in vehicle dynamics, there is a need to establish a coordinate
system. In this case, it was adopted an ISO 8855 coordinate system [9], where x is positive forward, y
is positive to the left and z is positive upward like shown in Fig. 2.1. It was used consistently throughout
this work to draw conclusions regarding forces and moments.
Forces are represented by F = (Fx , Fy , Fz ), being Fx in the longitudinal direction (negative mean
drag), Fy side force and Fz vertical force (downforce or lift when negative or positive respectively). Using
the same idea for moments, M = (Mx , My , Mz ) where Mx is rolling moment, My pitching moment
(positive nose down) and Mz yaw moment (positive turns the car to the left).
9
When appropriate, all results obtained are represented in dimensionless form. Not only does it
help to obtain/compare results at different speeds and air densities (temperatures) but it also eases the
process of image visualization. These dimensionless coefficients are also displayed multiplied by the
reference dimensions in a way that does not mislead conclusions and forces/moments can be easily
obtained when multiplying by the dynamic pressure q. For example CL .A and CMy .A.b where A and b
are reference area and length respectively.
2.1.2 Geometry
The geometry used for CFD simulations is a simplified CAD model of FST10e and it is derived from
the model used for aerodynamic simulations. This CAD is prepared to get a clean and simplified ge-
ometry for mesh generation, that include: adding/removing of all details that are not relevant to external
aerodynamics (bolts, lids, internal components, etc.), simplification of surfaces to reduce element count
where is not necessary (for example, uprights and rims), adding a driver and adjusting features to im-
prove cell quality. This model was made using Siemens NX® with the intention to parameterize and
speed up the interaction between CAD and CFD.
Suspension
The entire suspension was rebuilt to replicate the car motion depending on suspension geometry.
Suspension is moving about its pivot points, even though bellcranks are not physically represented, and
updates spring size based on the wheel travel. Regarding suspension parts, rims were further simplified,
uprights redesigned including brake discs simplification, it was decided to include suspension brackets
1
, it includes steering rack motion and other small adjustments were made to reduce cell count wherever
it was found appropriate. Special attention was given to the tires as they deform under load. When tires
are under vertical load, the sidewalls bend, lowering the suspension by decreasing tire loaded radius as
shown in Fig. 2.2 [10] .
Figure 2.2: How a tire carries a vertical load if properly inflated [10]
Tire shape, particularly closer to the ground, may change the flow structures created and affect its
wake. For example, an inconsistency between pressure taps (positioned in the diffuser) and wind tunnel
1 Becausethere is an acute angle between monocoque and some suspension arms, there needs to be a feature to ”help”
the prism layers to transit from one surface to another. Adding suspension brackets solves this issue (vs removing prism layers
in some section as previously) and gives a more accurate geometry to estimate wake downstream of suspension arms to the
underbody.
10
results was found on corner entry, result of the tire sidewall shape being different from one another in the
2012 Red Bull RB8 Formula One car. After improving the tire shape in CFD simulations, it revealed to
be a much higher source of low energy air to the diffuser, decreasing downforce as described by Adrian
Newey [11] . Bearing in mind this effect, some adjustments were made to the tire shape and position to
better replicate the tire.
To replicate this motion in CAD, it was implemented a variable called Loaded radius change for each
tire that lowers the wheel with respect to the ground. It adjusts the coordinate Z of the wheel suspension
points and updates tire geometry accordingly (Fig. 2.3(a)). Besides radius change, tire shape was also
a little tuned to reproduce bending of side walls. Also in Fig. 2.3(a), it is represented a perfectly round
tire in blue colour, while in red includes the adjustment made to the sidewall shape near the ground.
Lastly, it is also added a small plinth to each wheel at the contact patch to ensure a good mesh quality.
This contact region, creates an acute angle inducing meshing problems and poorly skewed cells. This
plinth is an extrusion to this intersected shape (Fig. 2.3(b) and its size can be changed in Expression
window of CAD geometry (Fig. 2.5(a)) .
Aerodynamic package
Regarding the aerodynamic package, trailing edges were thickened to increase the mesh quality of
prism layers. It was noticed that prism layers started to decrease in number when reaching the trailing
edge due to the airfoil thickness being close to zero at this point and the mesh could not adapt to this
angle (shown in Fig. 2.4(a)). To solve this issue, trailing edges in every wing were increased in thickness
by one millimetre to assist mesh generation while keeping the same chord. There was a very noticeable
difference with this modification, with the final results presented in Fig. 2.4(b).
The wings surfaces were split into trailing edge, lower and upper surfaces allowing specific mesh
controls to be implemented for each.
Adjustability
The entire CAD was constructed having in mind parametrization. Besides the five parameters in
study described in Sec. 2.1.3, several other features were included to be easily updated: static camber
11
(a) Trailing edge thickness = 0 mm (b) Trailing edge thickness = 1 mm
and toe angles for each wheel, loaded radius tire change, contact patch size, front wing height, rear wing
flap angles (DRS), suspension point coordinates and roll centre heights. Most of them are shown in Fig.
2.5(a). These features represent adjustments that can be made during car setup ( for example adding
shims to adjust camber angle) and others to easily update future prototypes if found fit.
(a) Siemens NX® Expressions window (b) Ride height change from 80mm to 20mm - Note the rear suspension arms motion
to update car geometry and difference in spring size
Ride height
As of many race cars, FS prototypes also take advantage of ground effect. For aircraft, this effect
reduces wing tip vortices generating less induced drag and also diminishes the upwash, but in race cars
as wings are inverted, ground presence will lower the pressure at the suction peak increasing downforce
generated. This effect is similar to that obtained for wings and can be explained with the method of
images [12]. The FST10e prototype is no exception, using this concept to generate more downforce
12
with less drag penalty. As ground proximity affects massively the behaviour of the flow around a wing,
the distance to the ground (ride height) is a key factor to consider in order to evaluate aerodynamic
performance.
Ride height is the distance between the lowest part of the sprung mass body 2 and the ground. The
static ride height can be changed by adding shims to the suspension push rod and dynamic ride heights
change due to the car motion. Accelerations and differences in speed will vary vertical loads, changing
the ride height along the track.
As each axle moves independently, this distance is commonly divided into front ride height (FRH),
and rear ride height (RRH) as displayed in Fig. 2.6. Both are measured at the centre of the car (y = 0)
and they define the vehicle’s pitch angle, which directly affects to the angle of attack of all wings,
RRH − F RH
θ = arcsin , (2.1)
W heelbase
where wheelbase is the distance between the centres of the front and rear wheels.
Roll
When a vehicle is cornering, it undergoes lateral acceleration that causes a change in its roll angle
(rotation about x-axis). This rotation results in different distances to the ground at each side of the car
and, once again, ground effect comes into play but in different ways depending on the side of the car. As
pressure fields change from side to side, both side force Fy and rolling moment Mx will be different from
zero. To take this behaviour into account, roll angle is included as one of the variables of study shown in
Fig. 2.7.
Roll axis is the axis around which the car rotates longitudinally and is determined by the suspension
geometry. As suspension geometry differs from front to rear, each has its own roll centre. Roll centre
(RC) is an imaginary point that can be found by intersecting the lines connecting the middle of the contact
patches to the instant centres. Instant centre (IC) is another projected imaginary point linking the two
suspension arms for each wheel. An example of construction of these 2 imaginary points is presented
in Fig. 2.8.
2 In a vehicle, sprung mass includes all the parts that are supported by the suspension.
13
Figure 2.7: Roll angle
The roll centre height (RCH) is an important aspect for vehicle dynamics because not only it dictates
the axis of rotation but also how much the vehicle rolls due to the existing arm between this point and
CG, producing a rolling moment. In fact, because roll centre depends on instant centre it is also an
instantaneous point for each suspension geometry. Although the term instant means at the precise
moment of time, it was found appropriate to simplify and assume the instant centres to remain fixed,
therefore roll centres and roll axis is maintained constant throughout the study.
The rolling axis is obtained by connecting the front and real roll centres, as shown in Fig. 2.9.
14
Steering
Formula style cars are open-wheeled and uncovered wheels are a major source of flow separation
3
regions and large size wake due to its considerable wake. This low energy air is produced not only
by the considerable wheel size and shape but also because of its rotation. There are flow structures
4
created by this phenomena [13] shown in Fig. 2.10, for example the tyre jet is one of the more
impactful vortices created for being close to the ground and entering the diffuser area. Although all of
these time dependant structures will not be captured in this work it was found appropriate to at least
mention the complexity around this region.
Figure 2.10: Iso-surface of Q criterion exhibiting coherent structures around a tyre wake [13]
The most efficient aerodynamic component of FST10e is the underbody (based on previous CFD
during the design process) that is positioned right downstream of the front wheels and it can reveal to be
quite sensitive to the wheel wake, which changes in shape, size and orientation based on the steering
angle. The steering angle is the angle between the wheel orientation and the longitudinal axis of the
car and they are not identical in both sides. Following an Ackermann steering geometry [1], inner wheel
steers more than the outer depending on the corner radius and, therefore, they are divided into two: δi
being inner wheel steering angle and δo outer wheel steering angle as represented in Fig. 2.11.
In addition to this, as the car turns differently when steering, the steering angle change is coupled
to the cornering radius R. When cornering, in the frame of reference of the car, incoming airflow is no
longer straight but rather curved (in circular shape assuming a steady-state condition). As free-stream
angle changes depending on position, it means side-slip angle β will also vary throughout the car. For
these two reasons - wheel wake and variable free-stream condition - steering angle is chosen as a
parameter of interest for this study.
To find each steering angle and its corresponding corner radius, data from a suspension kinematics
simulation provided by FST Lisboa vehicle dynamics department was used that consisted in rotating the
steering wheel one degree at each step, and then measuring steering angles and calculating the instant
turn radius.
3 This low energy air is commonly called dirty air in motorsport and it is associated with lower total pressure values, opposed to
15
Figure 2.11: Steering and yaw angles in cornering condition
Yaw
There are two factors that change the yaw angle: the side-slip angle and wind influence. Similarly
to an aircraft, β is the side-slip angle, measured between vehicle’s longitudinal axis and velocity (as
represented in Fig. 2.11). This change in direction is caused by the vehicle understeer and oversteer
behaviour.
Due to weather conditions, there is a presence of wind, that is related to the incoming airflow direction
and magnitude, differentiating car velocity relative to air from car ground velocity. As, once again, this
changes the free-stream velocity (magnitude and direction), it was considered relevant for this study and
yaw angle is another parameter of interest.
FS tracks layout is regulated as explained previously in Sec. 1.2.2, having all sort of turns and
straights being limited to a maximum of eighty meters. Cornering condition was chosen for a compre-
hensive analysis and numerical error estimation because it is the most frequent situation of a FS car and
to gather more information about this type of simulation that is used less often in FST. There are three
main elements of a CFD simulation workflow: pre-processing, solver and post-processing. The first two
will be explained in the following sections but the latter will be shown in deeper analysis in Sec. 3.2.
16
2.2.1 Pre-Processing
Domain
In a steady-state corner, the vehicle is turning with a constant radius. To model this situation in a CFD
simulation, it is imposed an angular velocity to the air and in order to align with free-stream orientation
the domain takes a shape of a partial hollow cylinder (or a donut sector) illustrated in Fig. 2.12. The
dimensions shown in this example were kept constant, resulting in different angles of revolution depend-
ing on the corner radius R in study. The domain extends about 3 times the length of the car upstream
and 10 times downstream. These values were chosen to minimise the influence of domain walls on the
results and are based on previous tests made within FST team.
Figure 2.12: Domain dimensions - Example for a 14.9 meter corner radius
Grid generation
In order to perform CFD simulations, it is necessary to get a discretized representation of the geom-
etry in study. To get a good surface quality, the geometry was cleaned up with a feature called Surface
Wrapper. Using this, it was possible to get a closed surface of the domain of interest (i.e. only the
volume between domain walls and the car) with all the gaps closed, non-intersecting surfaces and with
better quality due to the surface remesher associated with this feature. It makes sure the base geometry
for meshing is water tight, reducing the possibility of errors while generating the mesh. During this step,
interfaces are also created between radiator and fluid regions that is going to be further explained later
in Sec. 2.2.2.
After geometry clean-up, the next step was to generate the mesh for finite volume solver to be
applied. The type of mesh chosen to perform these simulations was the Polyhedral Mesher. In Star-
CCM+® , polyhedral mesh has been an optimization target, reducing mesh generation time in the latest
versions of the software. This type of mesh originates from a tetrahedral mesh, in which the centre of
tetrahedral cells are connected, giving rise to a polyhedral type mesh [14].
17
Polyhedral mesh are made out of polygon faces, increasing the number of surface cells, typically
around 14 cell surfaces for each volume cell. This feature increases the number of neighbour cells
(Fig. 2.13), allowing gradients to be more accurately calculated as they ”include” more information about
the neighbourhood and also reduces numerical diffusion effects due to having more faces in multiple
directions [15].
To refine areas of interest, some boxes were created to efficiently control the cell size. These areas
of interest are mainly to capture the car wake, and the developed flow structures (for example flow
separation, vortices generation and bursting). Starting with the wake regions shown in Fig. 2.14, three
different refinement areas were created:
1. The largest refinement region has a volume control size of around two hundred millimetres, extends
to the outlet and aims to capture the disturbed flow from the rear wing. Rear wing wake moves
upwards and this refinement area follows its height to the end of the domain.
2. Second largest region of refinement is coloured in blue and has a target volume size of around
eighty millimetres. Since the polyhedral mesh elements can increase in size rapidly, this box was
a way to smooth the transition from the more refined area to the coarser one.
3. The last of the largest boxes is represented in colour green having a twenty millimetres volume
control size. This box has a different shape to adjust to the rear wing wake without adding unnec-
essary elements upstream of it.
Closer to the car, two more areas of refinement were created with the aim of capturing the wake
generated by the front wing and downstream of it, mainly around the lateral diffuser:
4. Coloured in cyan blue in Fig. 2.15 and targeting around six millimetres, this box was created to
better resolve flow structures generated by the front wing. Vortices created in this area affect
downstream elements significantly but a grid to fully resolve them and avoid its diffusion is not
affordable.
5. The last box created for refinement was targeted for seven millimetres around the lateral diffuser
area. As this component is quite complex, it also generates vortices to create low pressure areas
and so, this box helps capturing all these structures.
18
Figure 2.14: Refinement boxes - In farfield region
At the wing trailing edges, reducing minimum size keeps its shape and allows prism layers to transi-
tion from top to bottom. A finer control was associated with the upper surface compared to the lower one,
to focus on the steeper adverse pressure gradients present in this surface without adding unnecessary
cells to the lower surface that can capture smoother favourable pressure gradients with fewer cells.
Mesh settings varied quite a lot throughout the components of the car. Twenty two custom controls
(excluding the five explained above) were created to regulate target and minimum sizes, surface curva-
ture refinement for leading edges, number of prism layers and respective height, and to disable prism
layers for some components. This customization was performed to decrease the computational cost of
each simulation keeping the numerical uncertainty at acceptable levels. This is the main advantage of
unstructured flow solvers and in particular of those that accept cells of arbitrary shape. In Fig. 2.16 it is
possible to see the volume mesh, and the different surface controls with different colours.
2.2.2 Solver
To model the physical phenomena involved in these simulations, mathematical models, boundary
conditions, material properties, numerical solution techniques and discretization schemes have to be
chosen that fit into the type of study. The options used in the CFD simulations in this work were:
1. Mathematical Models
19
Figure 2.16: Volume mesh and different named surfaces
• Steady: Time-averaging is applied to define the mean flow and so the flow is three-dimensional
but steady;
• Turbulent: Fully turbulent free-stream was chosen and this topic will be approached in Sec.
2.2.2 to check if there was a need to include laminar flow;
2. Boundary conditions
• All y+ Wall Treatment and Wall Distance: Shear-stress at the wall can be determined form
the definition or wall functions. More details are described in the following Wall Modelling
Section;
3. Material properties
• Air SSL: Incompressible fluid with density equal to ρ = 1.225 kg/m3 and a dynamic viscosity
of µ = 1.8.10−5 P a.s ;
4. Solution techniques
• Segregated Flow: Using this model, solver resolves the conservation equations of mass and
momentum sequentially, as opposed to coupled flow that solves them simultaneously and
converges faster;
5. Discretization techniques
Turbulence Models
First of all, an approach using RANS was the starting point. In a practical point of view, this is the
only approach that allows to perform simulations in a reasonable amount of time and it fits the flow
characteristics due to large region of separated flow and other strong viscous effects, having a low
Reynolds number (O(106 )) and approximately steady-state behaviour.
20
To obtain the RANS equations, it is used a Reynolds decomposition in Navier-Stokes (NS) equations,
that separates the flow variable into a mean and fluctuation,
φ = φ̄ + φ0 , (2.2)
where φ represents an instantaneous value (in this case for velocity components and pressure), φ̄ the
time averaged value, φ0 the flutuation around the mean value, and then an averaging of the equations
themselves.
Assuming a steady-state, incompressible fluid with no external forces applied, the RANS equations
can be stated as
∂ ūi ∂ ∂ ūi ∂ ūj
ρūj = −p̄δij + µ + − ρu0i u0j . (2.3)
∂xj ∂xj ∂xj ∂xi | {z }
Reynolds Stress
Although the equations went under an averaging process, it results in a term called Reynolds Stress
last term in RHS of Eq. (2.3) that includes averaged fluctuating velocity components. As fluctuations
are not directly calculated, this term has to be modelled using a turbulence model to get a closure of the
system of equations [18]. The turbulence model chosen was the k − ω SST, that is known for its ability
to handle adverse pressure-gradient flows when compared to other eddy viscosity models [16]. These
characteristics are suitable to external aerodynamics applications, in order to have a better separation
prediction.
The k − ω is a SST eddy-viscosity model with two modelled transport equations. One equation for
turbulent kinetic energy k and another for the specific dissipation rate ω. SST stands for Shear Stress
Transport, and it is responsible mitigate problems related to sensitivity of this model to free-stream
5
turbulence properties. It accomplishes this by changing to a k − model in the free-stream while
resolving near wall with the k − ω formulation.
Reference Frame
In the car’s frame of reference, while driving steadily in a constant radius turn, both air and ground
travels in a rotating motion. To mimic this effect, the whole domain was set into a rotating motion by
adding a rotating reference frame. This technique is commonly used in turbomachinery applications,
which enables the possibility of simulating constant motion of components with a steady-state model
without moving the mesh (stationary mesh). In this work, a reference frame is applied to the whole
region of study to simulate corner condition but it has other applications in car simulations such as to
model the rotating air between the rim spokes.
By adding this rotating frame, Navier-Stokes equations will be changed. The inertial velocity vector
~uI is defined by
~ × ~r ,
~uI = ~uR + Ω (2.4)
5ω = /k - relation between ω and turbulent dissipation rate
21
~ - Rotation vector; ~r - distance
where: ~u - velocity vector; I - Inertial notation; R - Rotating notation; Ω
vector from the cell to axis of rotation.
By applying this transformation to the incompressible, steady Navier-Stokes equations, in the rotating
frame and in terms of relative velocity (new unknown instead of absolute velocity) [19], yields
In addition to the convective, diffusion and pressure gradient terms from the original RANS equation
(Eq. 2.3), this transformation gives rise to two additional source terms: centrifugal and coriolis forces
(noted in Eq. 2.5). These terms are related to the inertial accelerations objects are subjected when
rotating. By adding these two source terms to every cell present in the region of study, a rotating domain
is simulated. The angular velocity of the domain varies depending on the corner radius in order to keep
the car tangential velocity and respective Reynolds number as constant as possible.
Wall Modelling
A good resolution of the near-wall region is essential to predict flow separation and boundary layer
growth over the car. The inner region of the boundary layer consists of three sub-layers that are repre-
sented in Fig. 2.17(a) where u+ is non-dimensional velocity and y + is a non-dimensional distance from
the wall to the first mesh node.
To capture all these sub layers accurately, prism layer meshing was included in relevant surfaces
of the car. As these type of cells increase substantially the cell count, the prism layer heights and the
number of layers were modified accordingly to the boundary layer growth. Target y+ (Fig. 2.17(b)) value
for these surfaces was around one, to be positioned inside viscous sublayer.
In some locations, it was found appropriate to disable the prism layers (coloured red in Fig. 2.17(b))
and use the traditional wall function approach, which approximates the boundary layer as shown in Fig.
2.18(b). This is called a high-y + treatment and allows an important reduction in near-wall cells. This
22
6
wall function approach was only used in suspension parts, wheels, main hoop that either have sharp
corners or separation zones should not vary significantly.
(a) Low y + approach- shear-stress calculated from the defini- (b) High y + approach - shear-stress calculated from the log-law
tion
The remaining surfaces of the car (aerodynamic package, monocoque and driver) had a more ac-
curate prediction of the boundary layer using a low-y + treatment (Fig. 2.18(a)) using a finer mesh.
Because both high and low y + were used, an All y + treatment was chosen to model the walls.
Transition Model
Boundary layer transition is the phenomenon of laminar flow turning into turbulent. Transition affects
boundary layer growth but also pressure fields due to laminar separation bubbles and reattachment that
can occur for the typical Reynolds number of the wings of a FS car. To evaluate the need of integrating
a transition model into the simulations, two simulations were compared, with and without a transition
model.
The model chosen was the Gamma transition model. This model is a one-equation model [21] that
adds a production term for the intermittency (zero for laminar flow and one for turbulent flow) [22].
It was observed that using the transition model added around two hours of computational time, while
forces, moments and the residuals did not show any significant change. A further analysis was made
to find the differences in the pressure field between both models. A static pressure coefficient delta
between both simulations is displayed in Fig. 2.19 that led to the conclusion that the biggest change
is mostly underneath the front wing, where it appears to have a laminar separation bubble followed by
reattachment leading to a reduction in pressure.
With these results, it was found the best option not to use the transition model as it increased com-
putational cost without a significant change in the relevant variables - integrated forces and moments
6A roll bar located alongside or just behind the driver’s torso to protect the driver during roll over.
23
Figure 2.19: Delta pressure coefficient to study the effect of transition
Cooling
In this section, it will be explained how the radiators and fans were modelled. As radiators increase
static pressure upstream and fans impose a pressure rise, these cooling parts were integrated into the
model, as they influence external aerodynamics. These components are positioned at the rear of the
car, as shown in Fig. 1.2.
To simulate the radiators, it was used the feature of Porous Region in Star-CCM+® . By setting a
porous viscous and inertial resistance, a radiator can be modelled by restricting the air into a specific
direction and imposing a pressure drop depending on the incoming airflow. By imposing these constants,
it allows to replicate the influence of a radiator without using a complex shape including all the fins,
reducing massively the cell count. To estimate these porous parameters, the fin geometry was simulated
at different velocities and by fitting a parabolic regression to the pressure drop versus velocity (work done
in FST) it is possible to obtain the constants [23]:
The fan positioned downstream of the radiator is modelled as a fan interface which imposes a pres-
sure rise depending on the local flow rate. This data is commonly displayed in the fan product datasheet,
and the data points from P-Q curve (Fig. 2.20) were introduced to simulate the respective pressure rise
for incoming air flow rate. During geometry repairing, it is also added a coordinate system in the middle
point of the fan to use as a centre of rotation for swirl.
Boundary conditions
The boundary and initial conditions used to model a steady-state cornering car were the following:
24
Figure 2.20: Fan curve (taken from NMB 12038 VA datasheet)
• Whole Domain: the entire region is set in motion with a rotating reference frame to simulate
steady-state cornering has explained previously;
• Inlet: A velocity inlet boundary condition is used to define free-stream condition. In this case it is
specified as V = 0 m/s with respect to the inertial frame so that the inlet has a constant angular
velocity visible in Fig. 2.21;
• Car: All surfaces of the car are no-slip walls in order to induce Vwall = 0 m/s to take into account
the real viscous flow behaviour, leading to a development of boundary layer;
• Wheels: The wheels have a wall rotation boundary condition. Each wheel has its own coordinate
system based on the wheel centre position to set the correct axis of rotation. The angular velocity
of each wheel is based on the wheel centre height, distance to the turning centre and steering
angle;
• Ground: The lower surface representing the ground is also a no-slip wall but this time is has a
velocity associated to simulate the moving ground. This velocity is the same as the inlet having a
constant angular velocity;
• Outlet: A pressure outlet boundary was used with a constant relative pressure of 0 P a to represent
the outflow condition. The value chosen is zero because it is relative to the atmospheric pressure;
• Walls: The remaining walls of the computational domain ceiling, inner and outer wall in which the
last 2 are hidden in Fig. 2.21 are modelled as slip walls. This results in the shear stress being zero
on these walls, to not develop unrealistic boundary layers;
In an attempt to reduce simulation time by allocating better the mesh refinements, it was used an
Adaptive Mesh Refinement (AMR) strategy instead of the boxes displayed in Section 2.2.1. Adaptive
Mesh Refinement is a method to adjust mesh size (refine or coarsen) based on a chosen criterion from
the solution. This means the refinement occurs during the simulation, updating the mesh depending on
the current solution. The criterion depends on the type of simulation and the flow characteristics that are
25
Figure 2.21: Boundary conditions - corner condition
meant to be better captured in the solution. In this case, the chosen mesh refinement criterion for each
volume cell is
pt − pref
S(r, θ, z) = | ∇ · ∇ 1 2
| × (AdaptionCellSize)2 , (2.6)
2 ρ(ωr)
| {z }
CpT
where pt is total pressure, pref is reference pressure, ω is the angular velocity of simulation, r is the cell
centroid radius (distance to rotation axis), CpT is total pressure coefficient and AdaptionCellSize is the
cell size in volume mesh that is two times the maximum distance between cell centroid and any of the
cell vertices [24].
In this study, the objective was to focus cell allocation to high values of the Laplacian of total pressure
coefficient. It was opted to use the divergence of a gradient instead of a gradient alone to prevent cells to
be allocated in areas with constant/high gradient values that can be easily modelled with a fewer number
of cells, especially with polyhedral type cells. This Laplacian value is then multiplied by the cell size by
recommendation of the software user-guide [24] to get an influence on the current size, to focus on the
bigger cells than needed instead of infinitely refining small cells. The power of two related to this cell
size was tweaked to find a good ratio between the Laplacian value and cell size.
The total pressure coefficient is calculated with a different dynamic pressure because the freestream
velocity is not constant throughout the domain. The better way found to fix this problem was to express
linear velocity as a function of angular velocity and corresponding radius with: V = ωr . Using this, it
was possible to achieve a total pressure coefficient of one in all free-stream. This dynamic pressure
adjustment was also used for post-processing plots.
After some experiments, the values found to better suit the case to refine the mesh are
26
refine if S > 1.5
Adaptive behavior keep if 0.025 ≤ S ≤ 1.5 . (2.7)
coarsen if S < 0.025
Moreover, it was specified a minimum adaption cell size of five millimetres, not to be applicable to
prism layers and an iteration frequency of one hundred - every one hundred iterations the solver runs,
an update is made to the mesh based on the criterion mentioned above.
A maximum refinement level of two was also set to limit the number of subdivisions of a cell. In a
polyhedral cell, subdivision can go up to 12-15 cells due to the high number of vertices as visible in Fig.
2.22.
After this setup, it was compared a simulation corresponding to a high-speed corner attitude with
AMR and the traditional refinement boxes where the results are shown in Tab. 2.2. The differences
between the two results are not drastic, being the pitch moment CMy the most significant change. It
was also possible to note that although the final mesh ended with a fewer number of cells, the total
computational time was increased due to the additional calculations and mesh division during this time.
As the time to generate the initial mesh was smaller, it took roughly the same total time.
In Fig. 2.23 it is presented a plane section view at a coordinate X = 475mm (intersecting the front
wing), displaying the total pressure coefficient and respective mesh. It is visible that the AMR approach
(Fig. 2.23(b)) refines cells to better capture the vortices structures, while in Fig. 2.23(a) the whole area
around the front wing is refined but only slightly, meaning that there a possibility some areas may not be
so well represented.
This same approach was then tried in a low-speed corner attitude, but it was not successful. Mesh
size increased drastically and fully converged simulation could not be achieved even when tweaking the
constants. As this criterion seemed too sensitive to the corner radius, it was found the best option to
27
(a) Refinement boxes (b) Adaptive mesh
Figure 2.23: Comparison of the effect of mesh refinement in the total pressure coefficient at x = 475 mm
discard AMR and continue with the traditional refinement approach for the remaining simulations.
In this section, numerical error derived from the models used to solve this problem will be broken
down into four distinct parts that will be estimated, to draw more founded conclusions when analysing
the results [25]: round-off error, discretization error, statistical sampling error and iterative error.
As these equations do not have an analytical solution, there is no possibility to estimate numerical
error by
φ = φi − φ0 , (2.8)
where φ is the error of a variable, φi is the approximate solution and φ0 the exact solution.
As computers have a finite precision, there is a possibility that the round-off error is significant but
to determine it exactly, infinite precision would be needed. To reduce this error, it was used a double
precision version of the software, i.e. 8 bytes of memory and has 14 digits of precision. Having this
number of digits, this source of error can be negligible when compared to the others.
Statistical error is usually linked to unsteady simulations. Even though this work is based in a steady
state assumption, oscillations of the results occur if there is a physical unsteadiness. In this case,
vortex shedding phenomenon is responsible for an oscillation in forces and moments. As expected, this
physical instability occurs most predominantly at the rear of the car, where vortices are being released
from the rear wing and diffuser. It was noticeable that, after a number of iterations, the results would
oscillate around a steady value periodically.
28
Standard deviation for the last 500 iterations of mesh convergence analysis is presented in Fig. 2.24
to exhibit the oscillation of around four Newtons in downforce caused by vortex shedding in the whole
car. Contrastingly, the standard error (or standard deviation of the mean) is negligible, that is hardly
visible in black. Simply averaging everything for a number of iterations allowed to reduce substantially
statistical sampling error.
Figure 2.24: Downforce standard deviation and error for the last 500 iterations
Initially, during mesh convergence study, 1700 iterations were ran and results averaged for the last
500. It was found that 1700 iterations were excessive, because the convergence was achieved much
sooner and then just oscillated around a steady value. For the remaining simulations, a total of 1000
iterations were ran, averaging all results for the last 250 iterations.
In addition to calculating the mean forces and moments, the same approach was used for field
functions, responsible for post-processing images. By comparing images from different simulations,
results could be misleading if comparing images in different phases of a vortex shedding cycle. A direct
side by side comparison is displayed in Fig. 2.25 with coloured total pressure coefficient and field lines
based on velocity vector using line integral convolution (LIC) 7 . While in the last iteration image some
vortices are visible downstream of the rear wing, they are no longer visible to the right presenting a
smoother field result of the averaging process.
Figure 2.25: Comparison of the effect of averaging in the total pressure coefficient
29
2.4.3 Iterative Error
An iterative method has to be used to resolve non-linear equations, inducing an error between the
exact and the approximate solution. Segregated solvers may increase this error as they solve equations
consecutively, but in this case a high number of iterations was ran in an attempt to minimize this contri-
bution. Using a coupled solver was not an option due to excessive RAM needed to resolve the equations
simultaneously. To monitor this error, it was also given attention to the residuals and variables of interest
convergence [27].
The discretization error is typically the greatest source of error and it derives from the fact that con-
tinuous functions are split into a finite number of volume elements to solve a set of governing partial
differential equations. Assuming this type of error to have a predominant effect over the others (which is
usually the case [27]), the numerical error can be estimated and assumed to be the same as discretiza-
tion error. A special attention was given to get a tangible value for the discretization error and get a
clear picture of what can be concluded. The remaining sources of error were seeked to be as small as
possible as will be shown next in this Section.
The discretization error can be estimated based on a mesh convergence study where results are
obtained for different values of mesh spacing parameter, h. This study was performed in a low speed
corner attitude (A1 - Attitude 1 / Skidpad) and then repeated for an attitude where only 1 parameter was
changed, visible in Tab. 2.3. In this case, the steering angle was changed and consequently cornering
radius, to represent a high-speed corner (A2 - Attitude 2) . This change allows to understand the trends
when varying a parameter in study.
The approach used to estimate the error was based on a power series expansion (following the
procedure in [28]). Using a typical cell size hi , convergence order p and constant α, the discretization
error φ can be estimated as
φ = αhi p . (2.9)
To obtain grid convergence, it is important that the grids are geometrically similar, meaning the
refinement should be approximately constant in all domain. To accomplish this in the best possible way,
a parameter was included in all relevant mesh controls (for example surface, volume sizes and prism
layers) to be easily changed and get the most similar refinement possible while using unstructured type
meshes.
30
Typical mesh size was varied in three distant intervals for each of the attitudes, obtaining the results
present in Tab. 2.4.
High-speed
corner
h = 1.0 37.9 31.6 1.629 0.084 -3.429 0.053 -0.212 -0.171
h = 1.4 15.2 14.0 1.626 0.083 -3.364 0.085 -0.196 -0.179
h = 2.0 6.9 6.7 1.588 0.078 -3.284 0.092 -0.183 -0.176
Besides presenting the results, it was found fit to also display the number of volume elements and
how long it took the solver to complete the simulation as they play an important role to make a decision
on which cell size fits best for this work. After obtaining these results, it was possible to estimate the
error and uncertainty. It turned out that the multiple variables convergence were non-monotonic, leading
to a calculation of two constants: α1 and α2 . Estimated uncertainty was proven to have a considerable
value, also due to the safety factor being three. These results can be seen in Fig. 2.26 for three variables
of interest downforce, drag and pitch moment.
For a clearer analysis, it is now presented Tab. 2.5, with respect to the lift coefficient to understand
the error and uncertainty orders of magnitude and the constants obtained in this process, where
φ0 or φ - Guess of exact solution or error for null cell size;
Uφ - Uncertainty;
α1 - First order error constant;
α2 - Second order error constant.
Having all this information, it is possible now to make a founded conclusion regarding the mesh size
to choose for further analysis. It was considered that results obtained using h = 1.0 did not justify the
31
(a) Downforce (b) Drag
Figure 2.26: Numerical uncertainty and error estimate for low and high-speed corner attitudes
extra computational cost when compared to h = 1.4. As this work is based on trends, practically the
same conclusions can be drawn while using h = 1.4 and saving half of the simulation time. Because a
high number of simulations is expected to be performed, it was found the best option for this work to use
h = 1.4 from this point onwards.
32
Chapter 3
Simulation Setup
In this Chapter, the remaining simulation conditions will be described, followed by the post-processing
data and finally a description of the macro created to automate the process.
Previously, in Chapter 2, cornering type simulation was analysed but more types of simulations were
performed in this work, namely straight line symmetrical and side wind simulations. These simulations
were particularly helpful for parametric studies where each parameter of study was isolated.
There is a case where symmetric condition can be assumed. When the turning radius (R), roll angle
(φ), yaw angle (ψ) and steering angles (δi and δo ) are null, it is expected that the solution is approximately
symmetric from one side of the vehicle to the other. Taking advantage of this condition, the domain of
study can be cut in half, reducing computational time by only studying half of the car and then ”mirroring”
the results to the opposite side.
To conduct symmetrical simulations, some adjustments were made with respect to boundary con-
ditions and domain of study in comparison to the corner type simulation, addressed in Sec. 2.2. The
domain is now shaped like a parallelepiped with the dimensions previously represented in Fig. 2.12
with exception of the width, that is now only six meters wide. In this type of simulation, there are no
additional source terms to the NS equations, and reference frames were not used since it is a straight
line condition. Only the right side of the car was simulated, as illustrated in Fig. 3.1.
The changes in boundary and initial conditions were:
• Whole Domain: An initial condition of constant velocity (equal to the inlet velocity) was applied to
the entire domain to speed up convergence;
• Inlet: A velocity inlet boundary condition was also used but, this time, set to V = 11 m/s in x-axis
direction, the air velocity in study;
33
Figure 3.1: Symmetrical simulation boundary conditions
• Wheels: The same approach was used to impose a wall rotation boundary condition but, this time,
it is only dependant on wheel centre height since there is no steering angle and the wheel speed
is constant;
• Ground: The lower surface representing the ground is also a no-slip wall, having a velocity asso-
ciated to simulate the moving ground. This velocity is the same as the inlet;
• Symmetry Plane: To model an imaginary symmetry plane, a symmetry plane boundary was
chosen. This type of boundary imposes a null value for normal gradients and variables of interest,
working the same way as a slip wall.
There are other attitudes corresponding to a straight line case (R = 0) that are not symmetric. These
conditions will be addressed in this Section.
In a straight line scenario, if either roll or yaw angles are different from zero, then a full car simulation
must be done. Although roll angle is not usually associated with straight line cases, it was useful to
include it here for parametric studies as will be discussed in Sec. 4.2. The yaw angle in a straight line
case is related to a presence of side wind, meaning that the air and ground velocities are different.
34
The setup of this condition is very similar to the symmetrical one. The boundary conditions that differ
are:
• Whole Domain: The domain is rotated around the CG so that freestream air orientation matches
the yaw angle in study ψ ;
• Inlet: The same procedure is applied here with V = 11 m/s, normal to the boundary surface to
keep the same Reynolds Number;
• Ground: Ground velocity is always applied in the vehicle’s x-axis direction, meaning that in the
side wind cases, the ground travels in different direction than the air, as shown in Fig. 3.2.
The remaining and last step of a CFD simulation workflow post-processing will be detailed, as item-
izing the images and results exported of each simulation and respective details. Lastly, a brief overview
of the macro created to automate CFD process will be given.
3.2.1 Post-Processing
The output files of each simulation will be detailed, namely: pre-checks, monitors, static and total
pressure coefficients. This last step is important to obtain the final results and to provide tools to analyse
and understand the mechanisms and changes that are causing those results.
Before heading into the specific data, a brief introduction about the views is now presented. In car
surface plots, different orientations were chosen for the purpose of acquiring all the information around
the vehicle. A total of fourteen views were created (all sides and corners), as shown in Fig. 3.3, and
they were used consistently throughout the entire work. In symmetrical cases, the results are mirrored
to the other side but, to not have repeated data, only nine images are produced (excluded Figs. 3.3(b),
3.3(d), 3.3(g), 3.3(k) and 3.3(m)).
(a) Bottom (b) Bottom Front (c) Bottom Front (d) Bottom Rear (e) Bottom Rear (f) Front (g) Left
Left Right Left Right
(h) Rear (i) Right (j) Top (k) Top Front Left (l) Top Front (m) Top Rear (n) Top Rear
Right Left Right
Figure 3.3: Views orientation of the surface plots around the vehicle
35
Pre-Checks
During the pre-processing phase, some images are produced for check-up. These pre-checks allow
to confirm/verify that the setup is correct and, if an error is identified, the simulation can be stopped
before the time-consuming process of solving the solution.
It starts with the geometry that, after being imported, it is plotted to check if every part was correctly
imported. An example of these images is shown in Fig. 3.3. Next, the coordinate systems are created
based on the surfaces imported (radiators, fans, wheels, CG) and different views are exported to verify
the auxiliary coordinate systems position, having one example shown in Fig. 3.4(a). After domain and
refinement boxes are created, the same procedure is followed, making images similar to Fig. 2.14. For
the latter stage of pre-processing, the surface mesh is verified to ensure that there are no parts inter-
secting the ground or any errors during mesh generation (Fig. 3.4(b)). The last step of this verification
process is to check the boundary conditions on the first iteration, shown in Fig. 3.4(c).
(a) Isometric view of coordinate systems (b) Left view of surface mesh (c) Top view of boundary condition
Monitors
During the solver process, it is important to track the convergence behaviour over the iterations using
several monitors. Residuals are one of the most important aspects to take into consideration to assess
convergence of a solution. It is a way to represent the imbalance resulted from the discretized equations
for all cells. Residuals were not scaled or normalized and the equations monitored were continuity,
momentum in each direction, energy, turbulent kinetic energy and specific dissipation rate. An example
of the residuals for a typical corner simulation is shown in Fig. 3.5 where they decrease some order
of magnitude and stagnate in around two hundred iterations exhibiting an oscillatory behaviour for the
solution.
Although residuals are fundamental to judge convergence, they are not sufficient. Other quantities
of interest must be also analysed such as forces and moments that were previously mentioned in Sec.
2.4.2. The quantities monitored during the solver convergence process were:
• Centre of Loads: This monitor was included to get an approximate value of the vehicle’s centre of
pressure. Centre of loads reports the position where the moment is minimal and it is a good way
to get a first idea of centre of pressure;
• Cooling: Pressure drop and rise in each radiator and fan respectively were monitored and also
mass flow rate across the cooling assembly to ensure the interfaces are working as desired;
36
Figure 3.5: Residuals monitor
• Forces: Forces in all three directions were monitored, including overall values and also breakdown
on each aerodynamic component. For the sake of easier visualization, two examples of force
monitors are displayed in Fig. 3.6, that only include the overall downforce and drag (excluding
force breakdown components);
• Moments: Overall moments were measured, in both CG position and origin;
• Solver Time: The last monitor was the solver time spent in each iteration. This helped in getting an
estimate time to complete the simulation and assessing if over usage of the computational power
was imminent.
Pressure Coefficient
The first post-processing images were related to surface (static) pressure coefficient Cp around the
vehicle
p − p∞
Cp = , (3.1)
q
where p is the static pressure and q represents the dynamic pressure that makes the coefficient adimen-
sional.
37
This allows to analyse the pressure field that ultimately gives rise to the resulting forces and moments,
better understanding the cause of the obtained results.
The dynamic pressure q defines the freestream conditions, or the undisturbed flow far upstream the
vehicle. There are two different ways to calculate this dynamic pressure depending on whether it is a
straight line using car velocity V , or a corner simulation where it depends on the local radius r,
1 ρV 2 ,
if R = 0
2
q= . (3.2)
1 ρ(ωr)2 , if R 6= 0
2
An introduction is now given to explain the legend structure created that is common to all images.
The image was divided to separate the image content at the top from the additional information about
the simulation at the bottom, as exemplified in Fig. 3.7.
As displayed in Fig. 3.7, the legend has four different regions of information:
• Colourbar: At the bottom left corner, it is shown the colour scale and bar used. Custom colourbars
were created for a better and faster analysis of the post-processing data. In the case of Cp , a
neutral colour (grey) is attributed to Cp = 0, to highlight areas where there is relevant change of
pressure. The pressure side is positioned to the right, and it is coloured by warm colour tones with
yellow and red. The upper limit has a value of one that is related to the stagnation areas and is
represented by a different colour (darker red/brown) to be directly identified. Low pressure side
has an increased range since the magnitude is higher, and is linked to cold colour tones. Since
green and blue were not enough to capture low pressure areas, pink was additionally introduced
to represent very low pressure areas;
• Attitude: Information about the simulation setup is displayed in two columns. The first column
includes front and rear ride heights, roll and yaw angles and the second column has the inner and
outer steering wheel angles, velocity magnitude and corner radius (or type of straight line simula-
tion). With this information displayed in every image, it is easier to understand the correspondent
attitude of each image when comparing side by side;
• Simulation ID: Each simulation had a unique run identification number. In the case of Fig. 3.7,
”MC” designate the initials of the author (Miguel Carreira) followed by ”074” that is the run number
of this work. Assigning an ID number to each simulation helps to organize and find the desired
data when needed;
38
• Axes: The last component of the legend is the display of axes to allow a better understanding of
the current orientation view.
An example of this image from the underside is now given in Fig. 3.8. Note the low pressure areas
in front wing in green and blue and even lower pressure in lateral diffuser in pink. It is also noticeable
the stagnation area in the front wing leading edge, monocoque’s nose and front wheels, represented in
dark red.
Figure 3.8: Example of static pressure coefficient Cp in bottom front left view
The second type of post-processing images is focused on displaying airflow information throughout
the vehicle rather than checking vehicle’s surface as previously. Cross sections are created for these
”slices” using constant coordinate planes and it is possible to see the flow structures by colouring total
pressure coefficient, related to the flow energy. It is possible to understand not only how vortices are
behaving during its formation or even bursting, but also separation zones and wake.
At first glance, these cross sections views might not be easy to locate so, some additional information
was attached to the legend. An example of a legend is shown in Fig. 3.9, where it is possible to see
two additional things that replace the axes information from Cp plots: plane and coordinate. In this case,
a side view of the vehicle is displayed in black, and with a red vertical line, the cross section plane is
represented. This red line updates depending on the plane coordinate that is shown to the right as in
the example ”X = −1565.0 mm” .
39
In Fig. 3.9, it is shown a section view legend in x-axis direction but the same method was applied for
the other two directions. The legend suffered minor changes for other cases: in y-axis, it is represented
a top view of the vehicle as in Fig. 3.10(a); regarding z direction sweep, the plane is now an horizontal
red line (shown in Fig. 3.10(b)).
Along the x-axis, an image is exported every fifteen millimetres, making a total of two hundred a
nineteen images for this sweep. Along y and z directions, the images are spaced by nine and ten
millimetres, respectively. With this information across the solution, it is possible to capture all the flow
features. An example of the information contained in each image is shown in Fig. 3.11.
As explained previously, the image is coloured by total pressure coefficient but there are also shown
subtle velocity field lines (LIC) to represent flow direction. A different colourbar was created for this
analysis to better adjust the flow characteristics. Undisturbed flow is represented in white and then
transits from warmer to colder tones, finishing in pink and then black to characterize regions of very
low energy / total pressure coefficient. This image is positioned between the wheels, sectioning the
expansion zone of the lateral diffuser. It is visible that on the left side of the figure the diffuser has
stalled, having a massive region with pink and black while on the right side airflow remains attached. It
is also noticeable the endplate vortices generated at the terminal plate of side elements, positioned on
both lateral ends, having the vortex core represented by black.
3.2.2 Macro
As a high number of simulations is expected to be made, it was found fit to automate the CFD
simulation to avoid a large amount of repetitive tasks, that are prompt to user error and consume time
40
that can be used more efficiently somewhere else. So, to better allocate time resources, it was invested
some time to create a macro that executes the entire CFD simulation workflow.
Simulation
With respect to the CFD simulation, the workflow was scripted having the goal in mind of just pressing
a play button to run the entire simulation and export the results for later analysis. The software used
is Simcenter Star-CCM+® version 2020.3 R8 (that uses double precision as mentioned before) and the
macros executed use the Java programming language. Several workstations were used to perform CFD
simulations, with 64Gb of RAM and a 4 or 12 core processor.
The first step to create the macro files consisted in performing each step manually while taking
advantage of the feature ”recording Macro”, that adds a command to the macro file for each action
made manually by the user. The files of each step were then compiled into several modules for each
major subtask, that are called by a more general macro that play the entire process. These files were
then adapted to integrate variables that are dependant on the input desired (the only interaction from the
user is to change the input variable in the code shown in Fig. 3.12(b)). The macro was also modified to
include all three types of simulation, based on the input by applying if conditions that are represented in
Fig. 3.12(a) and for cycles were useful mainly to export post-processing data. The macro workflow is
represented with a flowchart in Fig. 3.12(a), starting with the intent of creating a new simulation up to
exporting the results in a CSV format file and post-processing data.
41
Reports
As there is a substantial amount of data from all the simulations, an Excel sheet was made to gather
it all in one organized place for easy access when needed. To accomplish this, a small macro was pro-
duced in Excel that reads the exported CSV type file and copies the data into the correct cell depending
on the simulation. This macro was written in Visual Basic programming language and it was assigned
to the ”Import Data” button positioned on the top left corner of the datasheet, as shown in Fig. 3.13.
A section of the Excel sheet is shown in Fig. 3.13 that includes three example simulations to show
its structure. It starts with the simulation setup coloured in grey that includes the ID (project name and a
run number), a brief description of the simulation purpose, followed by the attitude in which the car was
set to. The simulation results are shown right after, with the forces, moments, efficiency and balance.
Following the overall results, there is a force breakdown in yellow containing both downforce and drag of
each component (some are omitted for brevity). In green, there is the cooling section, that has mass flow
and pressure drops and coloured in orange there is information regarding the solver and workstation,
such as computational time, number of cells, residuals, computer used and date of simulation.
Lastly, the moments at the origin, required to calculate the balance are shown. Two types of balance
are represented that consist in two moment balances about the wheel contact patch. Pitch moment
balance around front wheels (origin), My0 , gives rise to the front wheels ratio of downforce,
My0
Balance F ront [%] = 100 + × 100 , (3.3)
W heelBase × Downf orce
and a roll moment balance around the car centre, Mx0 , quantifies the proportion of load in the outer
wheels when compared to the total downforce,
Mx 0
Balance Outer [%] = × 100 . (3.4)
(T rack × Downf orce) + 0.5
42
Chapter 4
Parametric Analysis
Now that the setup and automation is set, all the tools necessary for study analysis are complete.
A parametric analysis is performed to evaluate the influence of each parameter of interest on the
results. It is also called sensitivity analysis, in which only one parameter is changed at a time. This study
includes all parameters of interest shown in Sec. 2.1.3: ride heights, roll, steering and yaw angles. A
baseline setup was chosen as a starting point for comparison as detailed in Tab. 4.1. Ride heights were
set at forty millimetres (front and rear), that is somewhere in between the ground clearance limits, using
a zero degree pitch angle as a default angle of attack. The remaining angles (yaw, roll and steer) were
also set to zero so that each variable can be isolated to minimize the influence on other parameters.
This attitude also allowed for the baseline to be in a symmetrical condition to decrease computational
time spent.
The first parameter in study is the ride height. As already mentioned, the front and rear ride heights
are connected to each other and can not be disassociated from one another. For instance, if we keep
the rear ride height constant and change the front one, the change in pitch angle leads to a different
expansion angle of the diffuser. For this reason, the front and rear ride heights will be analysed together.
The range of ride heights was set based on the travel suspension and usual ride heights configura-
tions, up to any part scraping the floor resulting in a range of [20 mm, 60 mm] for the front ride height
(FRH) and [20 mm, 80 mm] for the rear (RRH). As the rear suspension contains softer springs, the in-
terval of possible values is larger. A sample selection was needed to cover this entire envelope of ride
height values and a spacing of ten millimetres was chosen to evenly distribute the samples, totalling in
43
thirty-five simulations. All points correspond to symmetrical studies, speeding up this process.
The results obtained from these simulations are shown in Fig. 4.1, that display CL .A, CD .A and
pitch moment balance since the other forces and moments are zero, in these symmetric conditions. The
baseline (FRH40-RRH40) is coloured in white and then increase/decrease are represented in red or
green depending on the variable.
(a) CL .A (b) CD .A
Starting from the Fig. 4.1(a), an increment in CL .A is generally better so is coloured in green, while
red means a loss in downforce. There are two distinct areas, one around the bottom left corner in red
and at the top right corner in green. Lower values of front and rear ground clearance are associated
with lower downforce numbers. In an attempt to identify the cause of this behaviour, it was checked the
Cp plot for three configurations coinciding with the bottom left corner, image centre and top right corner
shown in Fig. 4.2 in the same order. It is visible that the low pressure areas at the front wing and lateral
diffuser are smaller in size and magnitude with lower ride heights. It was also noticed that the footplate
vortex [29] did not form when closer to the ground, leading to a loss in downforce in the front wing. The
low pressure associated with the vortex core is barely visible in Fig. 4.2(c) increasing low pressure area
under the front wing. With respect to the diffuser, it also reduces in efficiency when closer to the ground.
It might be ”saturated” and reached the point where downforce starts to decrease with lower distances to
the ground [30]. The less mass flow of air entering the diffusers may be also a reason for this behaviour.
The same overall trend can be found in the CD between bottom left and top right corners while
looking at Fig. 4.1(b), but this time it is displayed in opposite colours because lower drag is beneficial.
The induced drag is the reason these two variables are linked reaching the same conclusion.
However, the aero balance response was different. On the left side region in Fig. 4.1(c) the balance
shifts backwards, while on the right side has more downforce on the front wheels. It seems this balance
44
(a) FRH20-RRH20: CL .A = 2.805 (b) FRH40-RRH50: CL .A = 3.03 (c) FRH60-RRH80: CL .A = 3.411
Figure 4.2: Cp from bottom view for comparison between different values of ride heights
is more related to the vehicle’s pitch angle than the ground clearance itself.
The following parameter of interest is the roll angle φ, being the maximum value for this angle set to
two degrees. This value is not very high due to the high roll stiffness of this type of cars and together
with low ground clearance of diffusers and front wing. To cover this interval of [0◦ , 2◦ ], a step of half a
degree was chosen, keeping the remaining attitude settings at the baseline, giving rise to the results in
Tab. 4.2. For an easier understanding of the results trend, they are also represented graphically in Fig.
4.3.
It is possible to conclude that CD .A and CL .A have a downward trend with increasing roll angle.
As the outer side of the front wing gets closer to the ground, a similar behaviour to low ride heights is
observed, bursting the footplate vortex on the right side of the vehicle. This phenomenon is captured
45
Figure 4.3: Roll angle parameter study
in Fig. 4.4 where it is noticeable the burst vortex in Fig. 4.4(b), corresponding to φ = 1.0◦ , where the
wing is closer to the ground, but when φ = 0.0◦ in Fig. 4.4(a) the vortex is well generated. This bursting
reduces the low pressure region under the wing causing a significant drop in the front downforce when
reaching one degree (almost four percentage points), shown in black.
(a) φ = 0.0◦
(b) φ = 1.0◦
Figure 4.4: CpT section visualization for comparison between different roll angles
The balance starts slowly shifting backwards with increasing roll angle, being the only variable that
does not resemble a monotonic function. The reason for this change is that the right lateral diffuser
starts producing less downforce when closer to the ground. Mainly the flat outer portion of the diffuser
exhibits a smaller low pressure region, visible in Fig. 4.5. As this loss of downforce is located towards
the backside of the vehicle, it leads to an increase of the front wheel percentage load.
Since the left side is working ”better” than the right side, a side force CY .A is also produced in the
positive y-axis direction, as displayed with a green line. This increasing side force generates a negative
rolling moment (CM x ), transferring some load from the outer to the inner wheels, as seen in purple
colour.
46
(a) φ = 1.0◦ (b) φ = 2.0◦
Figure 4.5: Cp visualization from bottom front right view for different roll angles
The steering angle δ and consequent turning radius are the next variables targeted, with the averaged
steering angle varying within the range [0◦ , 12.5◦ ]. The upper limit value for this variable does not reach
the angle correspondent to the maximum steering wheel motion because of self intersecting domain
walls with low turning radius. With the full extent of the domain in a tight corner, the outlet wall would
intersect the inlet, making it impossible to simulate. Having this in mind, the highest value turned out to
be a bit less than seven meters of cornering radius. The sample points to simulate were spaced by two
and a half degrees of averaged steering angle between both wheels, since steering angles and corner
radius varied quite a lot. The outcome of these simulations are presented in Tab. 4.3 and then in graphic
mode in Fig. 4.6. Note the increasing difference between inner and outer wheel steering angle as it
increases due to Ackerman effect.
The drag and lift coefficient seem to remain fairly constant throughout the whole range of steer-
ing. There is a small increasing and then decreasing trend of the lift coefficient while drag is mostly
47
unchanged. The side force seems to increase linearly with the steer angle, in a substantial rise in mag-
nitude, starting from zero ramping up to three tenths. These changes in forces are due to differences
in low pressure regions located mainly at the front wing and lateral diffusers. The different sideslip an-
gles induced by the flow curvature across the vehicle alter the vortices behaviour and front wheel wake
change also interferes with the lateral diffuser performance.
The rolling moment CM x .A.b, expressed as aero balance of outer wheels (coloured in black) starts
to slightly increase at δavg. = 2.5◦ but then follows a declining trend, shifting the load to the inner wheels.
The yaw moment CM z shows a linear decreasing trend throughout the entire steering range. This
moment always increases in magnitude (with a negative sign) and, as z-axis points upwards, it induces
an understeer behaviour to the car. The flow curvature is the cause for this yawing moment, as the
different sideslip angles create a moment around the centre of gravity.
To illustrate the impact of the flow curvature onto the flow structures, Fig. 4.7 displays the evolution of
the iso-surface CpT = 0 across the steer angle interval, representing the wake throughout the vehicle .
The most visible differences are present in the rear wing endplate vortices and front wheel wake (in this
image right wheel is more visible). From 4.7(a) to 4.7(c), the rear wing right endplate vortex (left side of
image) starts to weaken and, in contrast, the left side vortex has strengthened together with a change in
trajectory, following the flow curvature. It is also possible to see that the right front wheel wake changes
in position and shape when steering, heading outwards. In the front wing right support plate 1 contains
an ever increasing in size vortex due to the side slip angle induced by the turning flow.
Figure 4.7: Iso surface of CpT = 0 evolution for different steering angles
The last variable of study is the yaw angle ψ. The interval is [0◦ , 15.0◦ ] and the points simulated
were separated by two and a half degrees of yaw angle, while maintaining the remainder settings of
the baseline attitude. Since considerable wind can occur during competitions, the limit chosen for this
variable was relatively high. The results are shown in Tab. 4.4, being displayed graphically in Fig. 4.8.
1 Support plate is the plate that serves as a mount from the front wing to the monocoque
48
Table 4.4: Yaw angle parametric study
Yaw angle Forces Moments
ψ[◦ ] CD .A CY .A CL .A CM x .A.b CM y .A.b CM z .A.b
0.0 1.510 0.000 -3.203 0.000 -0.358 0.000
2.5 1.497 -0.192 -3.119 0.057 -0.427 -0.007
5.0 1.498 -0.416 -3.060 0.052 -0.434 0.0004
7.5 1.447 -0.518 -2.930 0.084 -0.417 -0.041
10.0 1.384 -0.651 -2.828 0.086 -0.335 -0.067
12.5 1.347 -0.788 -2.693 0.130 -0.296 -0.093
15.0 1.317 -0.894 -2.674 0.125 -0.282 -0.132
This parameter of interest revealed to be one of the most sensitive ones, having a wide range of
results across the yaw angle interval. All force components feature a downward trend: the downforce
coefficient decreases in a somewhat linear way, having this time a much steeper reduction (almost
nineteen percentage points from the baseline), and the drag coefficient has the same behaviour, related
to less induced drag generated by less downforce. As expected, sideslip angle increases side force
substantially, reaching almost a nine-tenth increment of CY .A from the baseline. In the car’s reference
frame, the incoming air is coming from the left, so it makes sense that the side force has a negative y-
axis component. The front wheel load is the only parameter behaving in a non-monotonic way, starting
to drop down but, after ψ = 7.5◦ , it raises up. The aero balance shifts forward due to a sudden descent in
drag force on the rear wing. As the rear wing is positioned above and behind the CG, both the drag and
downforce contribute to a negative pitching moment, loading the rear wheels. As the drag decreased
significantly starting from ψ = 7.5◦ , it led to an increasing front wheel balance.
In Fig. 4.9, it is possible to understand the consequences of the yaw angle change into the pressure
side of the whole car. The driver head stagnation region gives an idea of the incoming airflow angle
and the rear wing high pressure reduced greatly at high yaw angles. It is also noticeable in Fig. 4.9(b)
some flat plates where the pressure is building up due to sideslip angle as: right side element endplate,
midplates 2 and support plates of front wing. The right side support plate even has a low pressure region
induced by separation.
Finally, the outer wheel load shows an ever increasing tendency, surpassing sixty percentage points,
shifting the centre of pressure outwards. Once again, as the incoming airflow comes from the left side,
the front left tyre wake enters into the diffuser, affecting its performance, while the front right tyre wake
is directed outward, meaning that the diffuser gets more energetic air comparing to the left side.
2 Midplate is positioned at approximately quarter span, to support the front wing flap and also reduce wing tip vortices
49
(a) ψ = 0.0◦ (b) ψ = 15.0◦
Figure 4.9: Cp visualization from top view from two different values of yaw angle
Overall, all the parameters studied presented a meaningful sensitivity either to some force component
or car balance, that is shown in Tab. 4.5 with percentage changes.
Even though the ride height only features three coefficient values different from zero (out of six),
they varied substantially throughout the entire working envelope. For example, CL .A varied within the
range [2.767, 3.548] reaching a delta value of 13 percentage points indicating a significant sensitivity to
this parameter. Drag and aero balance did not showed such great variation but they were also relevant,
changing eleven and fourteen percentage points, respectively.
The roll angle study has proved not to be so sensitive in magnitude compared to the ride height,
but was helpful to identify some monotonic trends regarding forces, outer wheel load and a peculiar
behaviour with respect to the front aero balance. The maximum magnitude change was the CL .A with a
nine percentage point decrease.
The third variable was the steering angle and a bit unexpectedly, it revealed minor changes in CL .A
and CD .A. Nonetheless, there were some other different important phenomena captured: substantial
increase in side force, aero balance shifting outwards by almost four percentage points and the major
change was the rise in yaw moment CM z .
Finally, the yaw angle exhibited a large interval of results, having the most impact of all variables in
drag and side force coefficient. CL .A consistently decreased within the interval [−3.203, −2.674] pre-
senting the largest delta of sixteen percentage points, and reached the maximum value for outside load.
Given all this data, all five variables were found relevant, so this work will continue to include them
all.
50
Chapter 5
Surrogate Model
This work aims to evaluate the aerodynamic performance of a prototype as a function of five vari-
ables. As the response of the design variables is found over time-consuming CFD simulations it was
opted to build an approximate model, also called surrogate model that attempts to reproduce the be-
haviour of the simulation model. This will enable a better understanding of the response over the gaps
left from the finite number of simulations in a more rapid and effective manner.
Constructing an approximate model requires several steps: design of experiments, construction of
the model itself and its validation. These will be approached in this Chapter, complemented by the
results obtained using this process.
Surrogate models are usually divided into two different applications: design optimization and de-
sign space approximation [31]. Surrogate model based optimization seeks to find an optimum value of
interest within the design space D, using an iterative search/update process, where additional simula-
tions/samples are ran around the predicted best value from the model to converge. This work, however,
focus on the other type of surrogate, design space approximation, in which the interest is only in the
overall behaviour of the system throughout the design space.
A surrogate model fb attempts to mimic the response of the ”exact” model f , focusing on understand-
ing the input-output behaviour, dismissing the proceedings needed to obtain the output (as a black box).
The symbol ˆ will be used whenever the variable in question is respective to the surrogate model. The
function f (x) depends on the design variables and the only way to collect data for it is through discrete
observations or samples. To cover the entire design space with discrete points, it requires choosing a
sampling plan x(1) , x(2) , . . . , x(n) that establishes an efficient spatial arrangement between samples
for evaluation through CFD simulations. This first stage is also called design of experiments (DOE) and
it focus on finding the input variables at each point in the input space [32]. In this case, x consists in five
design variables - front and rear ride heights, yaw, roll and steer angles. Each data point results in the
output variables y comprising of six different coefficients of interest for this work (Eq. 5.1).
51
CD .A
x F RH
1
CY .A
x2 RRH
−CL .A
x = x3 = ψ −→ y= (5.1)
Bal [%] − Out
x4 φ
Bal [%] − F ront
x5 δ
CM z .A.b
The main interest in choosing the points coordinates is to minimize the bias error. As the data is
acquired through a determinist CFD simulation, there is a numerical error associated with it. Besides
the bias error, the empirical error also includes a variance component that can be defined as [32] :
• Bias: quantifies the extent to which the surrogate model outputs fˆ(x) differ from the true values
f (x) calculated as an average over all possible data sets D;
• Variance: measures the extent to which the surrogate model fˆ(x) is sensitive to particular data set
D. Each data set D corresponds to a random sample of the function of interest.
Assuming that the data set includes some type of noise or random element for an average error with
a quadratic loss function, the expected error for bias and variance can be expressed as
n o2
Ebias2 (x) = EADS [fˆ(x)] − f (x) , (5.2)
h i2
Evar (x) = EADS fˆ(x) − EADS [fˆ(x)] , (5.3)
where EADS is the expected value considering all possible data sets.
These two quantities are linked to each other and there is a compromise between them. If a certain
model has low bias error, then the variance is usually high. It is possible to decrease both error com-
ponents by increasing the number of points in the data set but that adds computational cost. The bias
error can also be reduced by distributing the design points uniformly through the design space D [33].
Two methods that implement these strategies are Orthogonal Arrays (OA) [34] and Latin Hypercube
sampling (LHS) [35]. Although orthogonal arrays are able to have an uniform distribution, it has some
drawbacks: the OA may not exist under certain conditions and there is a possibility of including replicate
points, that would greatly increase bias error (bias error is usually more significant when dealing with
deterministic problems [33], as in this work). For these reasons, it was discarded and implemented
a Latin Hypercube sampling strategy that is also most frequent and can be easily implemented using
Matlab® function lhsdesign.
Latin Hypercube sampling is a stratified sampling approach with the restriction of each input variable
(xk ) having all portions of its distribution represented by input values [32]. Sampling is done by dividing
the range of each design variable xk into Ns (sample size) strata of equal marginal probability 1/Ns ,
(j)
and sample once from each stratum [35]. Denoting the sample by xk with j = 1, . . . , Ns and in this
case k = 1, . . . , 5 , the sample is then constructed with components of the various xk ’s matched at
52
random. This random element may prompt different levels of uniformity, that can be evaluated in distinct
ways: maximizing the minimum distance or minimizing the correlation between points for example. In
this case, it was found appropriate to use the correlation criterion that minimizes the sum of between-
column squared correlations. This way it is possible to focus on decreasing the dependence among the
variable inputs to get more information about the whole data spectre.
A total of fifty samples were distributed using LHS, adding to the existing simulations already made
for parametric studies, totalling in around one hundred observation points to construct the surrogate
model. This number of additional simulations was found to be a good compromise between simulation
time and dataset for the surrogate model. All sample points used are present in Appendix A.1, from
which starting on ID 62 up to ID 111 are the observations allocated using LHS.
There are two distinct types of surrogate models: parametric approach, that assumes there exists a
global functional form of connection between input and output that it is known; and the non-parametric
approach that utilizes different types of models according to the region, building up a global model (e.g.
Radial Basis Functions [36]).
This work adopted a parametric strategy, since the overall aspect is the main concern. The option
was to construct the surrogate based on a Polynomial Regression Model (PRG) due to its simplicity to
implement [37] . There are also other methods that would fit into this work such as Kriging Modeling
(KRG) [38] but it was decided not to enter in too much complexity at this stage, however it may be target
of a future and more intricate development in this area.
PRG is a technique that implements a regression analysis (hence its name) to obtain an estimated
value for the function of interest f , based on NPRG number of basis functions zj , with their respective
coefficients β. This approximation of the function f is essentially a Taylor series expansion. Intuitively,
this expansion is expected to have a better accuracy with increasing order of the polynomial m, therefore
different values for m were tested. For example, considering m = 1, this method expresses the linear
equation for each sample i and basis function j [32]
NX
PRG
(i)
fi (z) = βj zj + εi , E (εi ) = 0 , V (εi ) = σ 2 , (5.4)
j=1
where the error εi is assumed to be independent of the sample number with the expected value to be
null and variance being σ 2 . Representing the first part in matricial form we obtain
f = Xβ + ε , (5.5)
with X being a matrix containing all the basis functions for every observation point with a size of Ns ×
NP RG ,
53
(1) (1) (1)
1 z2 z3 ··· zNPRG
(2) (2) (2)
1 z2 z3 ··· zNPRG
X=
.. .. .. .. ..
,
(5.6)
. . . . .
(Ns ) (Ns ) (N )
s
1 z2 z3 ··· zNPRG
(i)
where the basis functions z1 are simplified with a number 1 as the first basis function is respective to
the constant term.
The basis functions zj include all combinations possible between the design variables coefficients
up to a determined maximum order m. For instance, a problem considering two design variables
(x1 and x2 ) fitting with polynomial up to order 2 gives result to the following basis functions: z (i) ∈
1, x1 , x2 , x1 x2 , x21 , x22 , meaning that there are NP RG = 6 coefficients to estimate.
The vector β̂ containing the unknown parameters can be estimated by
−1
β̂ = X T X XT f . (5.7)
In the end, after estimating the β̂ vector, the surrogate model approximated output can be found
simply by resolving
f̂ = X β̂ . (5.8)
The third step of making a surrogate model is validation, necessary to assess the approximate model
quality. This evaluation is possible by estimating generalization error, providing information to make a
decision about the best model choice.
There are several approaches possible to estimate the generalization error such as: Split Sample
(SS) and Cross-Validation (CV) [32]. The SS method is the simplest and, as the name indicates, the
data is split into a training and testing set. The training set is used to construct the surrogate model and
the testing set is only applied to compute the generalization error from the model. This approach does
not use all the computed data to build the model and the generalization error is highly dependant on the
dataset division. For these reasons, it was discarded and CV validation was chosen instead. There are
more complex ways to validate the model such as bootstrapping but it were not taken into account as
CV was found adequate for the time being.
Cross validation takes advantage of the same methodology of splitting the data, but enables the use
of the majority of the data (or even the whole data) as training set. The data is split into k-subsets
meaning that a surrogate is built k times, each excluding one of the subsets from the training data and
using it as testing, to validate the model. The surrogate generalization error is then estimated based on
the error from each of the subsets surrogate models.
The strategy chosen to select the data subsets was the leave-one out cross validation method, where
54
k is equal to the number of samples Ns , exemplified in Fig. 5.1. A surrogate is built Ns + 1 times, Ns
times to get all subsets into testing once, and one additional surrogate that uses all the data as training.
Each subset is taken into account as training data in Ns surrogate models and once as validation
data, being this one of the major advantages of the CV over the SS method. The main disadvantage of
this approach is that a high number of surrogate models has to be built to validate all the spectrum, but
since it is a computationally cheap problem to solve, it is resolved in a short period of time. In this case
Ns = 104 → (104 + 1) × 6 = 630 surrogate models are built for each polynomial order.
To estimate the generalization error using a mean square error (GMSE) we obtain
k 2
1 X (−i)
GMSE = fi − fˆi , (5.9)
k i=1
that represents the variation between the function of interest at the test point fi and the estimated value
(−i)
from the surrogate model that does not include i as training data fˆi . This generalization error was
calculated for each function of interest and for all orders of the polynomials to provide information to
make a founded conclusion on which order suits the problem best.
5.4 Results
The process of creating a surrogate model described previously was then executed for multiple poly-
nomial orders m. It was decided to perform the PRG for m ∈ N, m ∈ [1, 12] and computing the GMSE by
cross validation for each order to assess which order results in a smaller error. This value was chosen
based on Forrester et al. [40] saying that usually m ≤ 15, and a maximum value of 12 was found to be
plausible.
The results obtained from the regression analysis were consistent throughout all functions of interest
having the polynomial of order two the smallest GMSE in all cases as shown in Tab. 5.1. The error
reduces from order 1 to 2 and then it starts to increase almost exponentially up to the maximum order
chosen. Although the model is more ”flexible” with increasing order, it also enables the possibility of over
adjusting to the data noise and not capture the global solution. After a certain value for m (usually four),
matrix ranking started to be deficient and so, the estimated values for β̂ may not be the best. Together
with higher order terms, it may escalate the values, increasing substantially the error.
55
It should be noted that, for higher order polynomials, a considerably data set is required. which was
not available, and they might have led to poorly adjust models.
Table 5.1: Generalization mean square error for different orders of PRG
Max. order PRG NPRG GMSE
m β̂ −CL .A CD .A CY .A [%] Out [%] F ront CM z .A.b
1 6 0.0200 0.0013 0.0014 0.7194 7.3009 1.81e-04
2 21 0.0131 0.0006 0.0011 0.6979 6.0992 1.10e-04
3 56 0.0913 0.0021 0.0124 3.8691 30.1952 0.0012
4 126 10.696 0.1333 14.304 1.49e+03 1.86e+03 0.2220
5 252 16.117 1.0141 10.425 1.95e+03 1.87e+04 1.6458
.. .. .. .. .. .. .. ..
. . . . . . . .
12 6188 7.11e+09 1.70e+09 2.78e+06 2.39e+12 7.27e+11 5.11e+06
All results achieved are displayed in a matricial form in Appendix A.2 for all variables of interest. An
example of the obtained expression for −CL .A is represented in full extent by
ˆ
−CL.A = 2.0797
where the input units for ride heights is millimetres and angles is degrees. It is possible to see all the
parameters and respective basis functions with respect to the five design variables.
To understand how the surrogate fits into the sample data, the observations used for the yaw angle
parametric analysis were plotted against the estimation from the surrogate model. This example region
is represented in Fig. 5.2(a), where the points display the CFD output and the lines are the estimated
values from the surrogate model. It is possible to notice that the polynomial captures the meaningful
trends present in the sample data for this area of study. The results for the steering angle are shown
in Fig. 5.2(b)and it possible to see some variables do not fit so well, especially CL .A. The remaining
parametric results from the surrogate model are shown in App. A.3.
As a small portion of the whole design space is represented in the previous case, only displaying a
1D behaviour with constant ride height/roll/steer whereas the entire map is 5D. In an attempt to better
understand the global response throughout the domain, another type of plot was created. To accomplish
this, it was taken advantage of a chart tile plot (also called nested plot) where it is possible to include the
entire scope of four design variables. One of the design variables must be set as constant, in this case
56
(a) Yaw angle ψ (b) Steering angle δ
Figure 5.2: Comparison between results from surrogate and parametric studies
it was chosen to keep the yaw angle constant: ψ = 0.0◦ for Fig. 5.3(a) and ψ = 15.0◦ for Fig. 5.3(b) to
see both ends of the yaw interval. The only function of interest displayed is −CL .A as a demonstrative
example.
This plot is essentially composed by matrices of contour plots (tiles) in which each tile plot two
variables against each other. In this case, the row and column of a tile correspond to the front and rear
ride heights respectively (global vertical and horizontal axes). On each tile, the roll angle φ varies along
the horizontal axis and steering angle δ along the vertical axis that are represented on the bottom left
corner tile in Fig. 5.3(a) (F RH = 20, RRH = 20). Each of these two images include 117,000 points
calculated from the resulting polynomial expression Eq. (5.10).
In Fig. 5.3(a), it is possible to see some overall resemblance to Fig. 4.1(a), where warmer colours
(higher values) are positioned at the top right corner and cold (lower values) at the bottom left, meaning
the global ride height response seems to be similar. Regarding the other two variables, it is possible to
check the different behaviour of the roll and steer based on the ride height in study, both in magnitude
and shape. The bottom left peak has an almost circular shape while the top right, the peak moves to the
left and resembles more like an ellipse.
Looking at Fig. 5.3(b), it is noticeable that the response to ride height looks completely different than
before (note that the colourbar magnitude has changed). The model predicts higher values of downforce
at the bottom right area for a high value of yaw angle.
The CL .A evolution is displayed in more detail with more yaw angle intervals in App. A.4.
57
(a) ψ = 0.0◦
(b) ψ = 15.0◦
58
Chapter 6
On-track Validation
Aerodynamic validation is a rather unexplored field in FST Lisboa team in the past, and this part of
the work is aimed to start that process and try to find some correlation between on-track testing and
CFD results. My colleague Jaime Pacheco is currently correlating CFD to wind tunnel [41] and together
with this on-track testing, it promotes aerodynamic validation in both worlds for the future. As these tests
are unprecedented within the team and with the limited tools and conditions available, the focus will be
on trends rather than absolute values to estimate the relative difference between both worlds.
Throughout this Chapter, the on-track tests will be addressed, starting with the description of each
test, followed by the post-processing of the data gathered and finalized with the results obtained.
To estimate the car’s aerodynamic performance on track, two distinct tests were set: a constant
speed test executed to obtain the downforce and a coast down test to estimate the drag force.
The tests took place at Military and Technical Training Center of the Air Force, (”Centro de Formação
Militar e Técnica da Força Aérea, CFMTFA”), previously known as Ota Air Base [42]. Ota’s Air Base
which contains two runways (numbered in Fig. 6.1), where runway number 1 was the one used for these
aerodynamic validation tests. The runway with more than two kilometres long provides enough space to
drive the car steadily and perform the coast down from a high speed.
Figure 6.1: Military and technical training center of the air force, Google Maps
59
6.1.1 Constant Speed
The first test consists of running at a constant speed and estimating the downforce produced by
measuring the suspension deflection. There are several ways to measure the forces on the wheels,
including wheel force transducers and load cells in the suspension rods, but, with the material available,
it was decided to use the existing linear potentiometers located at each spring on FST10e to measure
the spring displacement due to the downforce generated. A steady-state is a necessary condition for this
test as any lateral or longitudinal acceleration also impact suspension travel ence the constant speed.
Knowing the springs stiffnesses and the corresponding displacement, it is possible to determine the
vertical load on each wheel by applying Hooke’s law. The current speed should also be measured to be
able to calculate CL .A.
Summarizing the constant speed test, the essential quantities to complete this test are total down-
force produced, air speed to calculate CL .A by applying Downf orce/q and dynamic ride height to get
information about the car’s attitude. The instrumentation used will be described in Sec. 6.1.3.
The most common way to estimate the drag force is by conducting a coast down test, which consists
in accelerating the car up to a certain speed and then coasting (stopping without using any power or
brakes) allowing the car to slow down rolling to a stop. The resistant force responsible for decelerating
the car is then determined from Newton’s second law where it is possible to isolate the part correspond-
ing to the aerodynamic drag.
Intending to express the vehicle motion in a simple equation, some simplifications must be assumed.
It is considered that the change in the loaded radius of the tyres is negligible and therefore the radius rd
is assumed to be constant throughout the coast down as there are no moments applied by the motor or
braking. The estimation of the parameter rw will be further described in Sec. 6.1.2.
where the force during coast down can be broken down into several resistance forces components.
Variable d represents the deceleration or negative acceleration during coast down, Rr , Rg and Ra are
the rolling, grade and the air resistance, respectively, and Rf is the resistance caused by the drivetrain
friction. The full deduction to this and the following equations is explained in full extent in [43] and only
some relevant and suitable steps are shown for brevity.
The mass map is the apparent mass of the vehicle, since there is an additional component to the
kinetic energy from the rotating parts of the wheels. Apparent mass can be expressed as
Iw
map = mt 1 + 2 , (6.2)
rw
as the car experiences a larger mass than the mass in translation motion mt , due to the extra component
associated with rotation relating moment of inertia Iw of every rotating part of the wheel and the loaded
60
radius rw .
Assuming a constant CD independent on the Reynolds number, the air resistance force related to
the drag coefficient, Ra by
1
Ra = ρCD A(v + w)2 , (6.3)
2
where, A represents the frontal area, v the car speed and w the wind speed. From Eq. (6.3), it is
possible to confirm that the drag force depends on the air speed and not solely on the vehicle speed,
so the wind should be taken into account for this study. With some manipulation and expanding the
remaining resistance and friction forces, it is possible to get a parabolic equation as a function of air
speed,
− a = C0 + C1 v + C2 v 2 , (6.4)
Figure 6.2: Least-squares approximation of experimental data (107 V-D points): d(v) = 0.135 +
0.00104v + 0.00025v 2 (CD = 0.342) [43]
Outlining the coast down test, the necessary measurements are air speed and car acceleration to
61
estimate drag coefficient from the parabolic regression term with CD .A = 2C2 map /ρ. Additionally, car
attitude during the test is once again required and will be addressed in Sec. 6.2.
6.1.3 Instrumentation
To measure all the required variables during the test, several sensors were used that will be briefly
outlined next.
Acceleration
To record acceleration during test, it was used the Attitude and Heading Reference System(AHRS),
whose purpose is to provide information about the three axes for multiple variables that include roll
pitch and yaw angles. In the FST10e, the AHRS is positioned close to the CG to measure the most
representative state of the car using the sensor MTi 670-series from Xsens [44]. This sensor, shown
in Fig. 6.3, includes a gyroscope, accelerometer and magnetometer (all 3-axis) but, in this work, only
the accelerometer was used to measure acceleration throughout the test. The standard full range of the
accelerometer is 10g that is more than enough for this application.
Velocity
The AHRS box containing this XSens sensor also includes an antenna with GPS but, unfortunately,
it was not operational at the time of testing to measure velocities at the same time. There was no
opportunity either to install a pitot tube to measure air speed that would assist in better accuracy. The
alternative used was to take advantage of the velocity estimator output from the control algorithm. There
is a device located inside the left diffuser that is responsible for the live data processing, called ETAS
ES910, shown in Fig. 6.4. With inputs from a variety of sensors placed all over the car, it is possible
to estimate all sorts of unknown variables by using an Extended Kalman Filter used by FST team [45]
that includes the longitudinal velocity of the car. Together with sensors, this estimated state is then used
to implement torque vectoring, traction control and yaw rate controller for example, with all this being
processed in ETAS. This estimated velocity was the one used to indicate the velocity of the car over the
test.
Potentiometers
A linear potentiometer was used for measuring the suspension positioning. Each wheel suspension
includes one of these sensors, positioned close to the spring and damper assembly, as shown in Fig.
62
Figure 6.4: ETAS ES910 [46]
6.5(b) for the rear suspension case. This sensor is part of the ELPM Series with 75 millimetre range
from the brand Variohm [47] which has motorsports application and shown more clearly in Fig. 6.5(a).
By analysing data from previous tests, it was noticed that the sampling frequency of these sensors
was around 10Hz, resulting in variations up to 3-4 millimetres from two consecutive measurements. As
these sensors are capable of a higher frequency of measurement, the sampling frequency was adjusted
to 1000Hz for a more accurate measurements. Such high frequency might induce some noise to the
data but it is easy to filter it afterwards if necessary.
(a) ELPM linear position (b) Installation in the rear left suspen-
sensor [47] sion
Additional Instrumentation
Some additional sensors were used, not directly to quantify something but to verify if everything was
working as expected. These sensors were located on the pedals and steering rack. Additional linear
potentiometers on the brake and throttle pedal allowed for confirmation if there was no applied force to
the pedals during coast down. The steering column contains a rotary sensor that gave information to
check if the car was not steered during runs and to identify the turn around between runs.
To make better use of a testing day as they are expensive, limited and very exhausting, the FST team
follows a procedure to provide orientation, organization and logistics of a testing plan.
Before testing, a plan was outlined for the testing day together with the responsible and additional
tasks were added to make use of the remaining time at the track including flow vis and wool tuft validation
for the aerodynamic department. It was also handled the logistics of the day by assigning a role to each
of the 12 people that would go to the test and others responsible for vehicle preparation before and
check it after the test.
63
During the day of testing, the vehicle and material for the day were already prepared but some
measurements had to be made before leaving the university. The battery and wings were mounted
and with the driver sitting in the car, the vehicle was put on top of four weight scales, one in each
corner/wheel, to measure the corner weights. Camber angle was also measured of all 4 wheels using
an inclinometer.
At the test site, all the material was unloaded and organized to be easily accessible, some tables and
a tent were mounted to check if everything is working as intended by live telemetry. After some pre-test
mechanical checks, cones were positioned along the runway to assist the driver by having some visual
references on where to accelerate and to start the coast down in both directions, to make sure it is
repeated at the same location. A few iterations were run to ensure there was enough distance to cover
roughly 10 seconds of constant speed and coast down in both ways. Power was limited to maximize the
battery usage as this test is not focused on accelerating rapidly and the motor speed was also limited
to 17.000 RPM. This speed is high enough to cause a considerable suspension displacement from
downforce but not excessive to run out of space during the steady time.
Right before runs were started, the ride heights were measured together with the tyre pressure on
all wheels. The first configuration was a medium distance of ride heights as the car usually runs.
The first two configurations of ride height were done in the morning, using all of the State of Charge
(SoC) of the battery. The last setup was then tested later in the afternoon, so the wind influence could
have been different associated with the different times of the day.
Having collected all the raw data during the test runs, the next step is to process it to estimate the
desired quantities.
Since there is no possibility to measure the ride height directly, an estimation route was taken to
calculate the dynamic ride height. It is relatively easy to convert the suspension travel to wheel travel
but the largest obstacle for ride height estimation is that the tyres deformation under vertical load can
not be neglected.
The FST10e tyres are Continental C19, which is a Formula Student type tyre with approximately 19”
diameter. The documentation of this tyre features the loaded radius behaviour as function of the vertical
load applied for 3 tyre pressures (IP-Inflation Pressure) and 3 different camber angles (IA - Inclination
angle), as shown in Fig. 6.6.
The tyre pressure during the test was 70kP a on all wheels which is not present in the documentation.
To solve this problem, it was taken advantage of the software WebPlotDigitizer [49], in which data points
were extracted from the Fig. 6.6 with approximately 70 data points for each of the 9 curves. From
this data, it was then found the best adjusting polynomial regression (which revealed to be of second
64
Figure 6.6: Tyre radius as function of vertical load [48]
degree) for each of the cases. Data was then interpolated to obtain the value corresponding to the
intended 70kP a pressure. The obtained functions to the testing pressure are visible in Fig. 6.7(a). To
estimate the loaded tyre radius during the test, two inputs are needed to the tyre model (Fig. 6.7(b)) that
interpolates between the curves for the camber angle and vertical load. .
(a) Calculated loaded radius vs vertical load for 70kP a tyre pressure (b) Tyre model diagram
From the potentiometer measurements, many dependant variables can be estimated through various
methods, as schematically shown in Fig. 6.8.
These are applied to each wheel individually and front and rear wheels are differentiated because
the suspension geometry is distinct. The first step is to convert the digital signal of the potentiometers
into its position between 0-75 millimetres. This allows to understand the suspension position at each
reading. Using a suspension kinematics simulation, also provided by FST Lisboa (similar to steering
one mentioned in Sec. 2.1.3) a wheel travel distance is imposed from -50 to +50 millimetres from the
default static position with a step of 0.1 millimetres, in which the behaviour of all suspension geometry
positioning is recorded. Although these results were integrated in the workflow with 1000 singular points
and used by a lookup table and interpolation between data points, it is shown graphically in App. B.5 for
65
Figure 6.8: Workflow for potentiometer data-processing
For the dependant variable potentiometer positioning, the independent variables include: wheel
travel, spring travel, change in camber angle and current motion ratio. Having found this, it is possi-
ble to compute the current wheel camber by adding camber change to the static value read before.
Spring force can be calculated by applying Hooke’s Law with the corresponding spring stiffness k and
spring displacement ∆xspring ,
with kf ront = 68[N/mm]; krear = 55.5[N/mm] . It is also defined a static condition that is determined
when the car is stopped right before the test that will be further explained in Sec. 6.3.1. Afterwards, the
delta values found are compared to the static ones.
The motion ratio M R quantifies the relation between the wheel and spring displacement,
W heel Displacement
MR = . (6.7)
Spring Displacement
Although motion ratio does not usually vary much through the entire working range of the suspension,
as this data was available it was also included and updated on each reading. As the displacements
are different, so are the forces and to find the respective force at the wheel, Fspring is multiplied by the
motion ratio.
To find the downforce, an additional component is also included due to lift generated by the wheel
assembly. As the wheel generates lift (approximately constant on all 4 wheels CL .Awheel ≈ 0.02) it has
to be discounted on the overall measured force at the wheel expressed by
1
Fwheel = Fspring .M R − CL .Awheel ρVair 2 . (6.8)
| {z 2 }
Wheel lift force
To estimate the last input needed for the tyre model, the vertical load can be found by adding the
66
corresponding corner weight that represents a large part of the vertical load,
Finally, all the information needed to estimate the dynamic ride height at each time is known and,
taking as example the front ride height, F RH is determined by the mean between left and right ride
height changes, since it is measured at the centre of the car,
where the subscripts F L and F R refer to the front left and front right wheel respectively, ∆x is the wheel
travel change and ∆r is the loaded radius change.
To change the static ride height during testing, shims were added and removed to the suspension
pushrod. The first configuration ran was the high ride height, followed by the low and then later in the
afternoon a medium configuration was defined, as shown in Fig. 6.9.
Figure 6.9: Front ride height configurations with different shim size
The static ride height is measured with a jig that essentially unscrews until it hits both the car and
ground at the same time. With this length set, it is then measured by a calliper to define the static ride
height. However, in the rear axle case, the distance measured is from the floor up to the monocoque ,
meaning this is not directly the distance defined as RRH in Sec. 2.1.3. By applying trigonometry to the
diffuser expansion ”triangle” it is possible to obtain the effective rear ride height by
tan(α) × d
RRHstatic = RRHmeasured − , (6.11)
cos(arctan( RRHmeasuredl −F RHstatic ) − arcsin( √ 2 tan(α)×d ))
l +(tan(α)×d)2
where α is the expansion angle (α = 3◦ , coincident with monocoque surface), l is the wheelbase (l =
1540) and d is the distance between the start of diffuser angle and rear axle (d = 741).
67
During the testing day, the measurement jig was lost and it was only possible to accurately measure
the last configuration with it. The remaining conditions were measured roughly to compare with the esti-
mated values to verify if they were plausible. Tab. 6.1 summarizes the different shim sizes used for each
configuration to calculate the change in pushrod length and to estimate the remaining ride heights. The
default suspension corresponds to F RH = RRH = 50mm, used in for the kinematic suspension stud-
ies. As the effective pushrod length changes with shims, this small change in the suspension geometry
was also taken into account.
To convert the shims size into ride height, it is necessary to understand the pushrod behaviour of the
suspension kinematics as previously only the spring to wheel travel was known. A simplified diagram
with the pushrod, bellcrank representation and spring is shown in Fig. 6.10, with the distances and
angles described in Tab. 6.2.
With this information, it is possible to calculate the relation between spring and pushrod displace-
ments (referred to as M Rspring−push ),
W × sin(θ)
M Rspring−push = , (6.12)
Q
The change in shim size is converted to a spring delta, which can be then translated into wheel travel
and subsequently to the missing ride height. This approach revealed to be successful as the estimated
68
values matched the rough measurements on track for the first two conditions. The three static ride height
configurations tested are displayed in Tab. 6.3.
The test took place on a sunny but also a bit windy day. It is important to determine the atmospheric
conditions because they influence the results obtained.
Air Density
Although there were no reading of the air density directly, it can be estimated by using ideal gas law
p = ρRT .
By looking at weather conditions around Lisbon it was possible to approximate the air temperature of
the day to a value around T = 24◦ C. Using Standard Sea Level (SSL) pressure of P = 101.325kP a and
with gas constant of R = 287.05J/kgK the air density for the test was found to be ρ = 1.1945kg/m3 .
Wind
The aerodynamic forces the car undergoes during the test do not depend on the car’s velocity but
rather on relative airspeed that is the sum of the two. As wind directly affects the airspeed, it is in the
best interest to measure or at least have an estimate of it. Ideally, a Pitot tube should be used to read
the air velocity but unfortunately, there was no possibility to include one. It was attempted to allocate
the testing to a low wind day, but the weather forecast is not reliable for long periods of time necessary
to plan everything and it revealed to be quite windy during testing day. We paid a visit to the weather
centre at the Air Base but the data collected by them could not be shared with us but revealed the wind
direction was relatively constant throughout the day, coming from the North meaning it was aligned with
the testing runway. Figure 6.11 shows a picture taken during the test with the car standing still before one
of the runs, where it is possible to see the wool tufts at the rear wing endplate all turned backwards due
to the considerable wind. The wind component revealed to be a massive adversity for the test results
and the greatest source of error for the on track testing.
The first approach to mitigate this problem was to average the runs both ways, as one had tailwind
and the other had headwind, to cancel the wind influence. When looking at the data, it was immediately
realized that the dynamic attitudes varied significantly from one way to another due to the increasing
aerodynamic load and this idea was discarded.
69
Figure 6.11: Wool tufts during testing day
The last approach used to try and solve wind speed unknown was to assume the same CL .A between
runs with similar attitude. Two configurations ended up with approximately the same dynamic ride height,
each in one direction. It was then assumed that the CL .A would be the same for the two runs, finding
the wind speed at which this would be verified. From this approach the wind speed is about 19[km/h]
used for the upcoming results.
Moving onto the data and results obtained from the on-track testing, firstly an overview of a run and
the data selection is explained. In Fig. 6.12 it is displayed an example of a round-trip for the first setup,
which contains 2 of the 6 runs performed for each configuration with the intervals used to obtain the
desired coefficients for each test. The complete data for some example runs is shown in App. B for
further analysis if desired.
On the left hand side, it is represented the reference static condition that establishes the car stopped,
in which variables are calculated relative to this position. The car is then accelerated up to the intended
speed and the first constant speed test run is executed for more than 10 seconds (confirmed on the
x-axis) where the data is averaged over this interval. Right after, it is performed the coast down up until
70
car speed reaches zero. Acceleration then goes down to a negative value as the car is decelerating.
The following step is to turn around the car, which is visible in the steering input. After some time,
the limit speed is reached again (17000 RPM from the motors) and the second run is executed, followed
by the coast down once again returning back to the begging where this process is repeated two more
times for each setup, totalling 6 runs.
The difference in downforce and car attitude between both consecutive runs is noticeable and signif-
icant and, for this reason, it was decided to isolate each direction and compare 6 different configurations
(from 3 setups). Together with the car speed, it is also shown the estimated air speed that is the mo-
tivation for this variation. In the test represented by number 1 there is a tailwind component in effect,
resulting in smaller airspeed than the car’s speed while in run number 2, air speed surpasses car speed
(that remains constant in both directions). One might find it strange that during coast down 1, air ex-
ceeds car speed at some point, but this derives from the fact air speed is shown an absolute value and
in that instance tailwind is larger than the car’s speed, meaning that in reality from the car’s reference
frame, the air is coming from behind.
Having implemented the plan of action explained previously by processing all the data, the results
from the track testing were found. The final experimental results of the constant speed tests are sum-
marized in Tab. B.1 and some data over run examples are displayed in Sec. B.4 to depict the differ-
ence/variance when repeating tests and comparing between different attitudes.
Previously, it was mentioned that two different runs were assumed to have the same downforce
coefficient to determine the wind speed. These runs are numbered as 3 and 4, in which the FRH has
changed 3.3 mm and RRH only 0.3 mm. This approach was found to be the best case scenario in this
work to obtain an estimative for the wind speed. It was assigned a number for each attitude, but the
numeration does not match the test sequence. This order was chosen to separate the wind direction
runs. From 1-3, the test has tailwind and 4-6 headwind. This mixing could mislead to wrong conclusions
coming from the simplistic way that the wind problem was solved, besides the fact wind is gusty inducing
even more error from this component.
The results are illustrated in Fig. 6.13(a) including uncertainty bars. As the wind speed inaccuracy
exceeds greatly all other sources of error, the experimental uncertainty was estimated by the change in
the downforce coefficient with and without a wind component, being approximately 1.0 change in CL .A.
There is a missing CFD data point in attitude 5 since the geometry intersected the ground. During this
run, the rear diffuser was a bit scrapped underneath confirming this theory that the car was too low.
Overall, what stands out is the significant difference between experimental and numerical results
across all attitudes. The gap is smaller in the first two attitudes but is still significant. The aim of this
aerodynamic validation was to evaluate the trends across different conditions of the car and not focus
on absolute values or differences between results. Having this in mind, the uncertainty magnitudes have
proven to be massive, meaning that even if some trend seems to be captured it may be misleading
71
(a) −CL .A results (b) Ride heights and respective standard deviation
because the trend may not be verified with other values within the standard deviation range. As all the
data is averaged across the time interval with constant speed, there are some variations over this time,
mainly related to the car attitude due to bumps or irregularities on the runway. To better showcase the
attitude change over this interval, Fig. 6.13(b) displays the averaged standard deviation of ride height for
each case. The red marker represents the averaged condition between the 3 runs (each one showcased
below in grey). The wind direction is also separated by marker shape (square vs triangle) enabling the
perception of the difference between consecutive runs with the same setup. From the figure, none of
the ride heights presents a significant variation, although RRH features higher values with a maximum
value of around 1.6 millimetres.
Analysing the first 3 attitudes (tailwind conditions), the decreasing trend and magnitude change be-
tween 1 and 2 looks similar in CFD and on-track results but it does not apply to attitude 3, diverging from
one other. Looking at the headwind conditions, a similar behaviour occurs, with the decreasing tendency
from 4 to 5 to resemble the track results (even though with a greater change) but then attitude 6 down-
force remains constant, while CFD results points to an increased value. Attitudes 3 and 6 correspond
to the runs executed at the end of the day after a few hours have passed since the other scenarios and
may have been submitted to different conditions (air temperature, wind gusts for example). This may be
the cause of the results deviation from the expected trend, unlike the remaining cases.
After processing all the data from the coast down it was possible to understand how the car behaved
during this test. An example of a coast down is shown in Fig. 6.14, which includes speed, downforce, ride
heights and longitudinal acceleration over the coast down time. From this figure, some perturbations are
noticeable causing some oscillations in the attitude, the most intense being around 511 and 512 second
time mark caused by some bumps on the track. Over this test, both ride heights vary substantially from
the beginning to the end, due to the loss of aerodynamic load compressing the suspension.
Taking the acceleration and plotting it against the air speed it is possible to obtain the quadratic
regression needed to find the drag coefficient. An example of this data is shown in Fig. 6.15, where the
72
Figure 6.14: Overview of coast down test data
raw data is shown against the parabolic curve with the best fit. This raw data is already filtered due to
the high level of noise coming from the accelerometer, visible in the last variable shown in Fig. 6.14.
The data from the coast down was treated in multiple ways to extract the best potential from it. It
was attempted to split the entire coast down range into 3 sections, and analyse each of them individually
as the attitudes vary considerably. This approach did not work as intended as it resulted in a greater
variation of the drag coefficient between sections, even reaching negative values of CD .A in some
particular cases. The other approach was to trim the data to exclude the first few seconds and the
last part of coast down. When the driver released the throttle, a transient state with the car pitching
nose down occurs, so the processed data was trimmed out to eliminate this influence. During the latest
section of deceleration, the linear and constant terms have a greater impact and, since the goal is to
estimate the second order term, it was found pertinent to exclude this latter part as it would induce more
errors in the drag coefficient (excluded values below 10 km/h) . This trim approach was found to be the
best suited for this case and it was taken for the coast down data.
By applying the same methodology to the 18 runs executed, resulting in the data described in Tab.
B.3 which indicates the averaged ride height in which the coast down was executed and the obtained
drag coefficient putting into practice the procedure explained before. Some additional measurements
and unknowns are shown in App. B.2, and a few example of the data throughout the run is displayed
in App. B.4. The data presents a high level of variation between the same setup runs. The third run of
configuration 6 was aborted due to electrical problems during coast down.
The experimental data is graphically presented in Fig. 6.16(a) including the results from the surrogate
model in the same attitude, uncertainty and standard deviation. The experimental bars represent the
73
CD .A standard deviation from the 3 sample points ran. The numerical uncertainty represented is the
one estimated in Sec. 2.4 regarding the mesh convergence. Interestingly, there is one experimental
point that falls below the surrogate one (attitude 3) which goes against the remaining trend. Usually, it is
expected to underestimate the drag due to all the simplifications made in the geometry and conditions
of the simulations. Once again, a major source of error may come from the air speed estimation just like
in the case of the constant speed test that greatly influences the outcome of these results. One should
note that conclusions must not be drawn between attitudes 1-3 and 4-6 as they represent opposite wind
directions. It is possible to see the high experimental standard deviation, due to the variation between
repeating runs. This means that the repeatability is poor, maybe due to track or weather conditions, for
example wind gust present in some runs. Other possible reason for this variation might be the fact that
the runway has a sudden inclination, having a longitudinal slope that goes against the assumption made
that the track is flat. It was also noticed through all the track data that the forces between right and left
wheels always differed from one another a bit. This might be an indicator that the track may have some
lateral slope as well, but no other cause was identified for this behaviour.
As stated previously, attitude changes over the course of coast down and Fig. 6.16(b) demonstrates
that variation. As the rear springs are softer, one would expect the RRH to change the most and that
is what happens. In the headwind cases, both ride heights variations are higher reaching a maximum
value of 4 millimetres for the standard deviation.
Again focusing on the trends instead of quantitative comparisons, it is found rather difficult to get
some certain conclusion as there is not a lot of confidence in the values obtained for this test and within
the uncertainty bars, trends could reveal to be the opposite as initially thought. It is possible to notice
that the change in CD .A is far greater in experimental results than in numerical. It can be seen that
surrogate values barely shift between attitudes. As the CFD results do not differ so much from the
surrogate model, it was found pertinent to not run additional simulations to draw the same conclusions
as with the surrogate model alone.
74
Chapter 7
Conclusions
The major objective of this work was to develop a deeper understanding of the aerodynamic be-
haviour of a race car, in particular for the FST10e prototype. A working envelope was defined and it was
studied the response of forces and moments coefficients in all three directions with respect to individual
and coupled operating parameters, representative of the car state. During this work, it was possible to
notice a change in more than 35% in downforce that confirmed the car behaviour dependence on the
attitude and condition proposition for this work.
Five parameters of interest were analysed that determine the car attitude. Roll angle revealed to have
a monotonic and rather independent behaviour resulting in 9% change in downforce and a shift of 12%
in pressure centre. Ride height individually presented a significant change in all response coefficients.
It varied by 10% in CD .A and 20% in CL .A between peak values in symmetric conditions, and it was
understood that contrary to what was expected, FST10e best working range corresponded to high ride
heights. The steering angle together with flow curvature, affects the most the yawing moment of the car
due to different incidence angles across the car. Finally, the yaw angle presented the highest sensitivity
to side force coefficient and also a considerable change in drag and downforce with 12% and 16% change
respectively. It was also noticed that the downforce coefficient was linked between yaw angle and ride
height, meaning that the optimal ride height zone changed based on yaw.
7.1 Achievements
1. It was possible to perform on-track aerodynamic validation for the first time by executing constant
speed and coast down tests. Even if not many absolute conclusions have been drawn it was a step
forward in the essential component of aerodynamic validation.
2. a further understanding of FST10e was gained, identifying the weaknesses and strengths of its
aerodynamic performance, helping in future design processes.
3. with an aerodynamic map, it was possible to quickly obtain results for any attitude by using a
surrogate model instead of spending many hours for CFD results.
75
To achieve the findings presented throughout this work, some tools were created:
1. A parametric CAD of FST10e prototype sped up the process of changing the car attitude and
provide a good surface quality for CFD simulations. In addition to the parameters of interest,
adjustability features were also incorporated to change the configuration if needed in the future;
2. a total of 116 CFD simulations were performed, focused on corner condition but also included two
straight line cases scenarios with a mesh convergence study to estimate numerical error;
3. Several macros were created to automate repetitive processes that are time consuming and prone
to user input error. The CFD workflow was automated and an Excel tool was made to allocate CFD
data into an organized document;
4. An algorithm to produce a polynomial regression surrogate model, that finds the best order based
on a cross validation leave-one out method;
5. A dynamic ride height estimator was developed to determine the car ride height during testing,
which also contains a tyre model to estimate its loaded radius.
Among the multiple subjects covered in this work, all of them present a large scope for future devel-
opment.
Numerical error and uncertainty from CFD results can be reduced to make more founded conclu-
sions. This can be achieved by increasing cell count or even by further exploring the AMR methodology
with other refinement criteria or even with other mesh type. This method was only briefly investigated
and it has a lot of potential to adapt the mesh to the problem in study. As computational power was
limited in this work and these simulations comprise of a complex geometry, they can be further refined
while running on faster hardware to decrease its error. Research can also be made by investigating other
turbulence models, for example resolving turbulent structures by the means of Detached Eddy Simula-
tion (DES) or even Large Eddy Simulation (LES) that are transient techniques in opposed to averaged
solutions present in this work. This approach can provide tools to study the time dependant structures
that caused the oscillations seen before but it also requires a high mesh resolution to be resolved.
The working range for the 5 parameters of interest can also be extended to incorporate a wider
envelope. This applies mainly to yaw angle ψ, where only positive values were taken into account.
These are related to oversteer scenarios but understeer conditions with negative yaw angles should be
considered in the future.
There are more complex approaches that can be used to construct the surrogate model and might
find a better response, with some examples being KRG and Radial Basis Functions. The surrogate
model can be integrated in future vehicle dynamics transient simulators within the team to update the
aerodynamic coefficients at each time step.
On-track testing validation is an area of focus that can be greatly exploited. With the possibility to
implement further instrumentation on the next prototype, correlation between CFD and reality can be
obtained using a quantitative approach estimating the absolute error from CFD.
76
Bibliography
[1] W. F. Milliken and D. L. Milliken. Race Car Vehicle Dynamics. SAE International, 1995. ISBN:978-
1560915263.
[2] J. Katz. Race Car Aerodynamics: Designing for Speed. Robert Bentley, 2nd edition, 2006.
ISBN:978-0837601427.
[7] X. Zhang, W. Toet, and J. Zerihan. Ground effect aerodynamics of race cars. Applied Mechanics
Reviews, 59(1):33–49, 01 2006. doi: 10.1115/1.2110263.
[9] I. O. for Standardization. Road vehicles — Vehicle dynamics and road-holding ability — Vocabulary.
International Organization for Standardization, Switzerland, ISO 8855:2011(E) edition, 2011. URL
https://www.iso.org/standard/51180.html.
[11] A. Newey. How to build a car. HarperCollins, 1st edition, 2017. ISBN:9780008196806.
[12] V. de Brederode. Aerodinâmica Incompressı́vel: Fundamentos. IST Press, 1st edition, 2014.
ISBN:978-989-8481-32-0.
[13] E. Croner, H. Bézard, C. Sicot, and G. Mothay. Aerodynamic characterization of the wake of an
isolated rolling wheel. International Journal of Heat and Fluid Flow, 43:233–243, 2013. ISSN
0142-727X. doi: https://doi.org/10.1016/j.ijheatfluidflow.2013.04.008.
77
[14] S. D. I. Software. Simcenter STAR-CCM+ User Guide, version 2020.3. In Polyhedral Mesher.
Siemens, 2021.
[15] F. Moukalled, L. Mangani, and M. Darwish. The Finite Volume Method in Computational Fluid
Dynamics: An Advanced Introduction with OpenFOAM® and Matlab®, volume 113. 10 2015. ISBN
978-3-319-16873-9. doi: 10.1007/978-3-319-16874-6.
[16] F. R. Menter. Two-equation eddy-viscosity turbulence models for engineering applications. AIAA
Journal, 32(8):1598–1605, 1994. doi: 10.2514/3.12149.
[17] V. Venkatakrishan. On the accuracy of limiters and convergence to steady state solutions. In 31st
Aerospace Sciences Meeting, 1993. doi: 10.2514/6.1993-880.
[18] G. Alfonsi. Reynolds-averaged Navier-Stokes equations for turbulence modeling. Applied Mechan-
ics Reviews, 62, 07 2009. doi: 10.1115/1.3124648.
[20] S. D. I. Software. Simcenter STAR-CCM+ User Guide, version 2020.3. In Wall Treatment. Siemens,
2021.
[21] F. Menter, P. Smirnov, T. Liu, and R. Avancha. A one-equation local correlation-based transition
model. Flow, Turbulence and Combustion, 95:1–37, 07 2015. doi: 10.1007/s10494-015-9622-4.
[22] S. D. I. Software. Simcenter STAR-CCM+ User Guide, version 2020.3. In Transition. Siemens,
2021.
[23] B. Çetin, K. G. Güler, and M. H. Aksel. Heat Exchangers - Design, Experiment and Simulation.
InTech, 2017. doi: http://dx.doi.org/10.5772/66281.
[24] S. D. I. Software. Simcenter STAR-CCM+ User Guide, version 2020.3. In Adaptive Mesh Refine-
ment. Siemens, 2021.
[25] W. L. Oberkampf and C. J. Roy. Verification and Validation in Scientific Computing. Cambridge
University Press, 2010. doi: 10.1017/CBO9780511760396.
[26] B. Cabral and L. C. Leedom. Imaging vector fields using line integral convolution. In Proceedings
of the 20th Annual Conference on Computer Graphics and Interactive Techniques, SIGGRAPH
’93, page 263–270, New York, NY, USA, 1993. Association for Computing Machinery. ISBN
0897916018. doi: 10.1145/166117.166151.
[27] L. Eça. Aerodinâmica Incompressı́vel: Exercı́cios. IST Press, 1st edition, 2015. ISBN:978-989-
8481-33-7.
[28] L. Eça and M. Hoekstra. A procedure for the estimation of the numerical uncertainty of cfd calcu-
lations based on grid refinement studies. Journal of Computational Physics, 262:104–130, 2014.
ISSN 0021-9991. doi: 10.1016/j.jcp.2014.01.006.
78
[29] J. M. Pegrum. Experimental study of the vortex system generated by a Formula 1 front wing.
Imperial College London, November 2006. doi: 10.25560/12585.
[30] J. Zerihan and X. Zhang. Aerodynamics of a single element wing in ground effect. Journal of
Aircraft, 37:1058–1064, 01 2000. doi: 10.2514/2.2711.
[31] A. C. Marta. Aircraft Optimal Design, chapter Surrogate Models and Multi-Fidelity. Instituto Superior
Técnico, Universidade de Lisboa, 2022.
[32] N. V. Queipo, R. T. Haftka, W. Shyy, T. Goel, R. Vaidyanathan, and P. Kevin Tucker. Surrogate-based
analysis and optimization. Progress in Aerospace Sciences, 41(1):1–28, 2005. ISSN 0376-0421.
doi: 10.1016/j.paerosci.2005.02.001.
[33] B. Tang. Orthogonal array-based latin hypercubes. Journal of the American Statistical Association,
88(424):1392–1397, 1993. doi: 10.1080/01621459.1993.10476423.
[34] A. Hedayat, N. Sloane, and J. Stufken. Orthogonal arrays: Theory and applications. Technometrics,
42, 11 2000. doi: 10.2307/1270971.
[35] M. Mckay, R. Beckman, and W. Conover. A comparison of three methods for selecting vales of input
variables in the analysis of output from a computer code. Technometrics, 21:239–245, 05 1979.
doi: 10.1080/00401706.1979.10489755.
[36] D. B. McDonald, W. J. Grantham, W. L. Tabor, and M. J. Murphy. Global and local optimization
using radial basis function response surface models. Applied Mathematical Modelling, 31(10):
2095–2110, 2007. ISSN 0307-904X. doi: 10.1016/j.apm.2006.08.008.
[37] N. Draper and H. Smith. Applied Regression Analysis. Wiley Series in Probability and Statistics.
Wiley, 1998. ISBN 9780471170822.
[38] J. Sacks, W. J. Welch, T. J. Mitchell, and H. P. Wynn. Design and analysis of computer experiments.
Statistical Science, 4(4):409–423, 1989. ISSN 08834237.
[39] N. Chlis. Comparison of Statistical Methods for Genomic Signature Extraction. PhD thesis, Techni-
cal University of Crete, 10 2013.
[40] A. Forrester, A. Sobester, and A. Keane. Engineering Design Via Surrogate Modelling: A Practical
Guide. John Wiley & Sons Ltd, 07 2008. ISBN 978-0-470-06068-1. doi: 10.1002/9780470770801.
[41] J. R. M. Pacheco. Wind tunnel testing of a complete formula student vehicle. Master’s thesis,
Instituto Superior Técnico, Universidade de Lisboa, 2022.
[43] I. Preda, D. Covaciu, and G. Ciolan. Coast down test – theoretical and experimental approach. In
CONAT 2010 - International Automotive Congress, 10 2010. doi: 10.13140/RG.2.1.4048.5925.
79
[44] MTi 600-series Data Sheet. Xsens Technologies B.V., Document MT1603P edition, 9 2019. Rev.
B.
[45] N. A. de Almeida Salgueiro. Torque vectoring control of an electric vehicle with in-wheel motors.
Master’s thesis, Instituto Superior Técnico, Universidade de Lisboa, 2021.
[47] Linear Position Sensor for Motorsports - ELPM Series. Variohm Eurosensor, 2021.
[48] Continental Formula Student Tire - Competition Tire 2019 (C19) Documentation. Continental, 4
2019.
80
Appendix A
x y
FRH RRH ψ φ δ Out Front
ID CD .A CY .A −CL .A CMz .A.b
[mm] [mm] [◦ ] [◦ ] [◦ ] [%] [%]
12 40.0 40.0 0.0 0.0 0.0 1.5097 0.000 3.2029 50.00 30.91 0.000
13 40.0 50.0 0.0 0.0 0.0 1.4604 0.000 3.0300 50.00 28.12 0.000
14 40.0 60.0 0.0 0.0 0.0 1.4670 0.000 3.2348 50.00 30.58 0.000
15 40.0 30.0 0.0 0.0 0.0 1.5013 0.000 3.2641 50.00 28.98 0.000
16 40.0 20.0 0.0 0.0 0.0 1.4571 0.000 3.0682 50.00 28.15 0.000
17 30.0 40.0 0.0 0.0 0.0 1.4478 0.000 2.9711 50.00 26.15 0.000
18 20.0 40.0 0.0 0.0 0.0 1.4377 0.000 3.0212 50.00 31.09 0.000
19 50.0 40.0 0.0 0.0 0.0 1.4825 0.000 3.1524 50.00 29.65 0.000
20 40.0 70.0 0.0 0.0 0.0 1.4807 0.000 3.2098 50.00 32.89 0.000
21 40.0 40.0 0.0 0.5 0.0 1.4988 0.003 3.1729 50.05 30.30 0.023
22 40.0 40.0 0.0 1.0 0.0 1.5011 0.018 3.0284 49.67 26.39 0.015
23 40.0 40.0 0.0 1.5 0.0 1.4735 0.038 3.0074 49.32 26.80 0.020
24 40.0 40.0 0.0 2.0 0.0 1.4536 0.051 2.9145 49.16 26.98 0.032
25 40.0 80.0 0.0 0.0 0.0 1.5290 0.000 3.3558 50.00 35.15 0.000
26 40.0 40.0 2.5 0.0 0.0 1.4974 -0.192 3.1190 53.06 29.10 -0.007
27 40.0 40.0 5.0 0.0 0.0 1.4978 -0.416 3.0598 54.85 28.58 0.000
28 40.0 40.0 7.5 0.0 0.0 1.4470 -0.518 2.9296 56.83 28.48 -0.041
29 40.0 40.0 10.0 0.0 0.0 1.3837 -0.651 2.8280 58.32 30.10 -0.067
30 40.0 40.0 12.5 0.0 0.0 1.3472 -0.788 2.6932 61.38 30.45 -0.093
31 40.0 40.0 15.0 0.0 0.0 1.3167 -0.894 2.6737 62.31 30.89 -0.132
32 40.0 40.0 0.0 0.0 2.538 1.4936 0.065 3.1978 50.51 29.16 -0.102
33 40.0 40.0 0.0 0.0 5.083 1.5097 0.120 3.2286 49.48 33.06 -0.201
34 40.0 40.0 0.0 0.0 7.6445 1.5175 0.190 3.2465 48.47 29.84 -0.295
Continued on next page
81
Table A.1 – Continued from previous page
x y
FRH RRH ψ φ δ Out Front
ID CD .A CY .A −CL .A CMz .A.b
[mm] [mm] [◦ ] [◦ ] [◦ ] [%] [%]
35 40.0 40.0 0.0 0.0 10.232 1.5195 0.275 3.2311 47.87 29.73 -0.404
36 40.0 40.0 0.0 0.0 12.854 1.5081 0.334 3.1740 46.67 30.22 -0.487
37 30.0 30.0 0.0 0.0 0.0 1.4632 0.000 2.9878 50.00 24.90 0.000
38 50.0 50.0 0.0 0.0 0.0 1.4974 0.000 3.5480 50.00 35.94 0.000
39 60.0 60.0 0.0 0.0 0.0 1.5450 0.000 3.5195 50.00 33.85 0.000
40 20.0 70.0 0.0 0.0 0.0 1.4565 0.000 3.0786 50.00 32.37 0.000
41 50.0 20.0 0.0 0.0 0.0 1.5014 0.000 3.1765 50.00 25.03 0.000
42 30.0 20.0 0.0 0.0 0.0 1.4921 0.000 3.0432 50.00 28.24 0.000
43 20.0 20.0 0.0 0.0 0.0 1.4126 0.000 2.8047 50.00 22.79 0.000
44 20.0 30.0 0.0 0.0 0.0 1.3950 0.000 2.7666 50.00 25.91 0.000
45 20.0 50.0 0.0 0.0 0.0 1.4454 0.000 3.0169 50.00 31.50 0.000
46 20.0 60.0 0.0 0.0 0.0 1.4871 0.000 3.1048 50.00 32.90 0.000
47 30.0 50.0 0.0 0.0 0.0 1.4627 0.000 3.0235 50.00 29.54 0.000
48 30.0 60.0 0.0 0.0 0.0 1.4499 0.000 3.3319 50.00 35.83 0.000
49 30.0 70.0 0.0 0.0 0.0 1.4546 0.000 3.2647 50.00 36.60 0.000
50 50.0 30.0 0.0 0.0 0.0 1.4969 0.000 3.2363 50.00 27.46 0.000
51 50.0 60.0 0.0 0.0 0.0 1.5452 0.000 3.4237 50.00 34.72 0.000
52 50.0 70.0 0.0 0.0 0.0 1.5378 0.000 3.3999 50.00 33.02 0.000
53 60.0 20.0 0.0 0.0 0.0 1.5328 0.000 3.2457 50.00 23.52 0.000
54 60.0 30.0 0.0 0.0 0.0 1.5451 0.000 3.2451 50.00 24.19 0.000
55 60.0 40.0 0.0 0.0 0.0 1.5661 0.000 3.4254 50.00 28.83 0.000
56 60.0 50.0 0.0 0.0 0.0 1.5604 0.000 3.3264 50.00 29.57 0.000
57 60.0 70.0 0.0 0.0 0.0 1.5628 0.000 3.4957 50.00 35.62 0.000
58 20.0 80.0 0.0 0.0 0.0 1.5007 0.000 3.0193 50.00 29.28 0.000
59 30.0 80.0 0.0 0.0 0.0 1.4852 0.000 3.3101 50.00 37.29 0.000
60 50.0 80.0 0.0 0.0 0.0 1.5402 0.000 3.3291 50.00 33.09 0.000
61 60.0 80.0 0.0 0.0 0.0 1.5532 0.000 3.4112 50.00 35.99 0.000
62 46.8 37.4 3.75 0.9 9.626 1.4856 -0.010 3.2795 49.60 29.38 -0.355
63 54.0 59.0 12.15 0.1 5.372 1.3376 -0.670 2.7305 59.46 31.41 -0.260
64 30.0 71.0 7.35 1.5 10.6215 1.4588 -0.201 2.9615 51.47 27.79 -0.404
65 55.6 32.6 9.15 1.6 7.875 1.3698 -0.417 2.8119 56.60 20.74 -0.286
66 29.2 24.2 3.15 0.98 5.1965 1.4634 -0.070 2.9774 51.72 27.76 -0.194
67 24.4 50.6 10.65 1.34 12.127 1.3733 -0.406 2.6675 54.73 27.66 -0.487
68 36.4 72.2 6.15 0.82 11.1265 1.4788 -0.199 3.2163 53.73 38.55 -0.406
69 38.8 48.2 11.85 1.86 8.127 1.3641 -0.538 2.7802 57.69 30.51 -0.332
70 59.6 56.6 11.25 0.94 9.3715 1.3502 -0.505 2.7072 57.01 27.37 -0.376
71 51.6 23.0 11.55 0.38 8.626 1.2898 -0.552 2.9494 55.77 24.04 -0.355
72 32.4 74.6 8.25 1.46 1.8775 1.4119 -0.448 2.9408 55.73 31.65 -0.083
73 58.0 33.8 0.75 0.7 4.373 1.5295 0.041 3.2720 48.00 27.40 -0.149
74 56.4 47.0 4.35 0.5 7.124 1.5392 -0.133 3.2566 50.76 28.21 -0.283
75 38.0 61.4 14.25 0.74 2.8755 1.3399 -0.834 2.8217 62.19 38.31 -0.195
76 52.4 65.0 3.45 0.46 10.8745 1.5144 0.002 3.2568 49.98 33.25 -0.401
77 27.6 77.0 10.95 0.02 1.6255 1.4170 -0.643 3.1789 57.46 40.48 -0.152
Continued on next page
82
Table A.1 – Continued from previous page
x y
FRH RRH ψ φ δ Out Front
ID CD .A CY .A −CL .A CMz .A.b
[mm] [mm] [◦ ] [◦ ] [◦ ] [%] [%]
78 34.8 57.8 8.85 1.78 2.3715 1.4156 -0.493 2.7627 56.71 29.79 -0.100
79 54.8 35.0 14.85 1.94 1.374 1.2491 -0.848 2.2898 64.13 23.02 -0.121
80 42.8 31.4 1.05 1.02 2.6235 1.5144 -0.023 3.4335 50.34 31.69 -0.088
81 53.2 30.2 7.05 0.18 6.1225 1.4059 -0.380 3.1443 53.44 26.69 -0.237
82 41.2 68.6 1.35 0.3 4.6225 1.5696 0.042 3.6285 51.53 38.18 -0.187
83 22.0 26.6 14.55 0.78 11.622 1.2384 -0.619 2.3787 58.35 25.47 -0.490
84 58.8 79.4 0.45 1.82 6.3765 1.5835 0.215 3.3616 49.45 36.07 -0.234
85 50.0 25.4 10.35 0.06 1.1225 1.3073 -0.672 2.7055 58.76 26.55 -0.096
86 46.0 67.4 12.75 0.26 11.374 1.3686 -0.560 2.8539 57.04 33.13 -0.482
87 44.4 38.6 13.95 0.54 12.379 1.2483 -0.623 2.6030 58.15 26.24 -0.509
88 40.4 36.2 1.95 1.42 6.621 1.5212 0.053 3.3126 49.54 31.72 -0.236
89 28.4 41.0 6.45 1.18 0.625 1.4217 -0.455 2.8452 54.98 29.54 -0.024
90 42.0 62.6 5.25 1.22 3.375 1.5181 -0.266 3.2016 55.56 33.55 -0.139
91 25.2 44.6 12.45 0.86 5.625 1.2948 -0.652 2.7111 58.78 33.42 -0.274
92 45.2 43.4 8.55 1.54 8.878 1.4147 -0.361 3.0094 55.65 28.36 -0.327
93 22.8 63.8 13.35 0.34 3.622 1.3525 -0.733 2.8513 58.49 35.76 -0.223
94 31.6 60.2 4.65 1.62 9.126 1.4801 -0.054 3.0770 49.73 29.72 -0.314
95 21.2 49.4 5.55 0.42 9.8715 1.4561 -0.164 3.1876 52.22 35.11 -0.373
96 26.8 78.2 2.85 0.6 2.125 1.5216 -0.114 3.2658 50.95 32.83 -0.080
97 43.6 39.8 1.65 0.22 0.875 1.5120 -0.114 3.2489 51.24 31.47 -0.030
98 48.4 66.2 13.65 1.26 4.1275 1.3208 -0.760 2.7176 61.87 35.76 -0.208
99 33.2 69.8 9.45 1.3 4.875 1.4213 -0.478 2.9218 55.95 33.53 -0.214
100 57.2 54.2 10.05 1.66 5.872 1.3997 -0.499 2.8250 57.98 28.54 -0.239
101 20.4 45.8 2.55 0.14 7.625 1.4674 0.023 3.1728 49.05 34.45 -0.307
102 30.8 29.0 13.05 1.7 0.375 1.2907 -0.753 2.6679 60.45 31.31 -0.072
103 23.6 21.8 9.75 1.3 3.8755 1.3510 -0.521 2.6766 58.88 27.08 -0.156
104 49.2 53.0 6.75 1.14 10.375 1.4371 -0.210 3.0319 53.18 29.48 -0.392
105 50.8 73.4 7.95 0.66 0.125 1.5081 -0.532 3.1240 59.04 35.14 -0.046
106 34.0 75.8 7.65 1.06 11.871 1.5115 -0.256 3.1179 53.60 34.16 -0.448
107 26.0 20.6 0.15 0.62 10.125 1.4720 0.248 2.9699 47.12 27.89 -0.350
108 35.6 55.4 4.95 1.1 7.375 1.5409 -0.157 3.3469 53.02 34.55 -0.267
109 37.2 51.8 2.25 1.38 8.375 1.5528 0.105 3.4450 49.97 35.90 -0.319
110 39.6 27.8 4.05 1.58 3.125 1.4734 -0.233 3.1720 52.38 30.28 -0.089
111 47.6 42.2 5.85 1.74 6.8755 1.4828 -0.199 3.1793 53.51 28.41 -0.254
112 38.3 22.0 1.0 0.39 0.0 1.4737 0.075 3.1661 49.53 27.32 -0.011
113 45.0 35.4 1.0 0.326 0.0 1.5086 0.070 3.2118 48.57 29.64 -0.005
114 49.9 43.0 1.0 0.46 0.0 1.5410 0.089 3.4017 49.68 31.80 -0.007
115 48.3 35.0 1.0 0.783 0.0 1.5194 -0.102 3.1828 50.62 28.10 -0.026
116 43.3 27.3 1.0 0.655 0.0 1.5055 -0.101 3.2036 50.47 27.61 -0.014
83
A.2
T
1.3225 0.0092525 2.0797 50.3087 18.0373 0.019904 1
0.0040983 0.00040121 0.031531 0.03227 0.25895 −0.00068258 F RH
Results
0.0009701
−0.00045281 0.0085836 −0.036386 0.15816 −0.00035564
RRH
0.0034308 −0.074534 0.026217 0.82502 0.29307 −0.0039085 ψ
0.018779 −0.001836 0.1034 −1.6456 1.6198 0.033951 φ
0.0066618
0.024799 0.048496 −0.054187 0.926 −0.038518
δ
−7 −5
−4.716 × 10−6 −4.5532 × 10 −8.1955 × 10 0.0010873 0.0024456 1.6048 × 10−6
F RH × RRH
−5
CD .A −0.00020473 −9.4447 × 10 −0.00093097 0.0033418 −0.01263 6.504 × 10−5 F RH × ψ
0.00063677 0.00023588 0.0019702 0.035123 0.060052 −0.00056083
F RH × φ
CY .A
−5
−0.00015694 1.1562 × 10 −0.00026495 −0.0046363 −0.0091112 3.0529 × 10−5
F RH × δ
−CL .A
−5 −5 −5
=
7.2088 × 10 −1.3983 × 10 −2.3425 × 10 −0.001248 0.0020487 −9.4857 × 10−6 × RRH × ψ
84
Bal [%] − Out
−0.00036556 0.00054645 −0.0029023 −0.0032019 −0.060223 −1.9296 × 10−5
RRH × φ
Bal [%] − F ront −5 −5
8.0103 × 10 −1.9935 × 10 0.00020373 0.0018843 0.0076417 −2.2927 × 10−5
RRH × δ
CM z .A.b
−0.00074869 −0.0004287 −0.0025126 0.08357 −0.035523 0.00062124
ψ × φ
| {z }
f̂ −0.00059635 −0.00042217 −0.00084686 −0.0063881 −0.032995 0.00031132 ψ×δ
0.0021642
0.00051953 0.0069771 −0.011353 −0.011244 −3.9408 × 10−5
φ × δ
2
−1.5815 × 10−5
−4.5226 × 10−6 −0.00022188 −0.0010983 −0.0045911 6.9853 × 10 −6
F RH
−6 −6 −6
−6.6173 × 10−7 3.7769 × 10 −3.9984 × 10 −9.7591 × 10 −0.00087591 3.025 × 10−6 RRH 2
−0.00079359 2
0.0012339 −0.0018687 −0.0048986 0.018499 −0.00043352 ψ
2
−0.019863
−0.0029673 −0.089235 −0.13373 −1.5823 0.0024393
φ
−5
−0.00017183 4.1506 × 10 −0.0037686 −0.0095627 −0.069302 4.0068 × 10−5 δ2
| {z } | {z }
β̂ X
(A.1)
A.3 Parametric Results
(a) Front ride height - RRH = 40mm (b) Rear ride height - F RH = 40mm
Figure A.1: Comparison between results from surrogate and parametric studies
85
Figure A.3: Yaw angle ψ = 2.5◦
86
Figure A.5: Yaw angle ψ = 7.5◦
87
Figure A.7: Yaw angle ψ = 12.5◦
88
Appendix B
Validation Results
Run −CL .A
[mm] [mm] [N ] [N ] [N ] [%] [km/h]
1 49.89 43.01 643.36 163.78 479.58 25.46 2.803 70.57
2 49.79 42.56 670.22 170.22 500.00 25.40 2.922 70.55
1
89
B.2 Additional Information
Moment of inertia of each wheel package Iwheel = 148243101.7[g.mm2 ] measured in CAD include
tyre, rim, brake disk, centre lock, hub and transmission parts.
Run CD .A
[mm] [mm] Fr-RH Rr-RH
1 52.11 52.87 1.418 1.00 1.85
2 52.07 52.78 1.995 1.14 2.06
1
3 Aborted
Avg. 45.77 41.26 1.616 1.67 3.82
90
91
(a) Setup 1 Run 1 (b) Setup 1 Run 3 (c) Setup 1 Run 5
Figure B.1: Full runs for Setup 1 - Tailwind runs corresponding to attitude 1
92
(a) Setup 1 Run 4 (b) Setup 2 Run 4 (c) Setup 3 Run 4
Figure B.2: Example Runs 4 corresponding to headwind condition
93
(a) Setup 1 Run 3 (b) Setup 2 Run 3 (c) Setup 3 Run 3
Figure B.7: Suspension kinematics based on wheel travel for ±50mm heave
97
98