An investigation on the switching of asymmetric wake flow and the bi-stable flow states of a simplified heavy vehicle

This study investigates the effect of aspect ratios (Ra*) on the wake bi-stability behind the GTS model using improved delayed-detached eddy simulation (IDDES) at a Reynolds number of 2.7×104. Eight cases with Ra* in the range of Ra*∈[0.5, 2.0] have been considered to investigate the asymmetric wake flow and the corresponding aerodynamic response. The primary objective of this study is to reveal the mechanism of the Ra* effects on the GTS’s flow topology, asymmetric wake flow and the bi-stable flow switching process. The accuracy of the numerical method has been validated by comparing the recirculation bubble, vortex core and aerodynamic drag with those obtained from the previous large-eddy simulation study and the water channel experiment. The results show that three typical flow topologies (FT-I, FT-II and FT-III) are classified according to the drag coefficient and the length relation between the vertical and horizontal recirculation region. The vertical and horizontal asymmetric characteristics can be transformed between two adjacent flow topologies. Moreover, a switching process of the bi-stable wake flow is also observed as Ra* is equal to 0.71. Additionally, the ground motion can reverse the dominated flow state and inhibiting the corresponding switching process.


Introduction
In energy conservation and emission reduction, the aerodynamic drag has become one of the most notable aspects in transportation. During the progress of the aerodynamic drag reduction of bluff bodies, the asymmetric wake flow structures are observed behind symmetric square-back models, such as the Ahmed body (Ahmed et al., 1984;Grandemang et al., 2012;Östh et al., 2014) and the ground transportations system (GTS) (Ghias et al., 2008;Salim & Ong, 2013;Unaune et al., 2005), which thereby has attracted a lot of attention from scholars in the field of vehicle aerodynamics. Furthermore, previous studies found two possible flow state solutions in the near wake of the GTS model and Ahmed body, where the two possible flow state solutions are called as the bi-stability flow.
Previous studies have found that the asymmetric characteristics of the turbulent wake flow behind the GTS model in the vertical direction are drastically impacted by the numerical method and mesh (Rao et al., 2018a;Rao et al., 2018b) as well as the gap between the fixed floor and the GTS's lower surface in the water channel tests CONTACT JiabinWang wangjiabin@csu.edu.cn (McArthur et al., 2016). Moreover, the effect of aspect ratios on the asymmetry of wake flow in the vertical direction has also been reported in the study of Hassaan et al. (2020) and Zhang et al. (2022b), and the results show that the aspect ratio is an essential geometric factor that affected the turbulent wake flow. Nevertheless, most of the previous studies concerning the GTS's aspect ratio effect are limited to the model with the value of Ra * larger than 1.0 and wake topology classification according to the asymmetry characteristics in the vertical plane. The incomplete exploration of the Ra * 's influence on the wake structure, especially for the Ra * value below 1.0, will certainly contribute to an uncomprehensive classification of the GTS's wake topology. Since the wake flow topology is closely related to the aerodynamic drag coefficients of the vehicles (Zhang et al., 2022b), one of the significant purposes of this paper is to establish a complete research system for the GTS's wake flow topology, to provide a more comprehensive understanding of the developing trend of GTS's aerodynamic drag, which thereby guides the vehicle engineers on designing new configurations for a square-back vehicle.
The effects of Ra * on the asymmetric characteristics of flow structures in the horizontal plane behind the Ahmed body and its corresponding aerodynamic behaviors have been widely investigated (Mcqueen et al., 2014;Venning et al., 2015). For the square-back Ahmed body, the aspect ratio value is equal to Ra * = 0.74, the asymmetric characteristic of the wake structures is extremely dominant in the horizontal plane, and the wake's asymmetric characteristic is observed in the vertical plane when the value of Ra * increases to 1/0.74 = 1.34 (crucial aspect ratio value of vertical asymmetry behind the squareback Ahmed body) from Ra * = 0.74 (crucial aspect ratio value of horizontal asymmetry behind the square-back Ahmed body) (Grandemang et al., 2013). This indicates that the occurrence of the wake's asymmetric characteristic of the square-back bodies might transform between vertical and horizontal planes in case the aspect ratio value changes to its reciprocal. Besides, current research concerning the Ra * effects on the wake flow asymmetry mainly focus on the case of Ra * ≥ 1.0. For the investigation of GTS's aspect ratio value lower than 1.0, the Ra * value will certainly reach the reciprocal of the crucial aspect ratio value of vertical asymmetry the GTS model reported by McArthur et al. (2016). While, at the current stage, the appearance of the wake's asymmetric characteristic of the GTS model in the horizontal plane for the Ra * value lower than 1.0 has not been reported yet.
The switching process of the bi-stable asymmetric wake flow state has been extensively studied for the turbulent wake flow behind the Ahmed body. The random switching of the bi-state flow states can be found in long-time wind tunnel tests (Volpe, Devinant, & Kourta, 2015) and even relative short-time numerical simulations (Dalla et al., 2019;He et al., 2021). Few studies observed the switching of the asymmetric wake flow state. The existing study on the switching process of the bi-stable asymmetric wake flow state behind the GTS model is limited to the vertical middle central plane, which is induced by the numerical scheme. As the GTS's Ra * decrease close to Ra * = 0.74 [the aspect ratio value of the square-back Ahmed body reported by Grandemang et al. (2013)], the switching of the asymmetric wake flow state might occur in the horizontal plane. Nevertheless, to the authors' knowledge, the corresponding switching process of the bi-stable asymmetric wake flow state in the horizontal plane has not been reported yet.
In summary, the objectives in this study are summarized as follows.
• To perform a comprehensive investigation of the Ra * effects on the GTS's turbulent wake structure and corresponding aerodynamic performance with a wide Ra * range from 0.5-2.0, which thereby provides guidance to the vehicle engineers when designing new configurations for a square-back freight vehicle. • To identify whether the wake's asymmetric characteristic of the GTS model will transform between the vertical and horizontal planes, therefore providing a more comprehensive cognition of the GTS's wake asymmetry. • To investigate the switching of the asymmetric wake flow state in the horizontal plane and to reveal ground motion effects on this random switching process, thereby providing theoretical guidance for the experimental investigation of the bi-stable flow in the wind tunnel.
The remainder of this paper is organized as follows. Section 2 provides the numerical method setup and validation. In section 3, the Ra * effect on the GTS's aerodynamic drag, flow structure and asymmetric characteristics of wake flow and the ground motion influence on the random switching of the wake flow state in the horizontal plane have been analyzed. The conclusions are detailed in Section 4.

Numerical set-up
The original GTS model is presented in Figure 1. In this paper, the width W (W = 0.064 m) of the GTS model was defined as the characteristic length, and the length (L) and height (H) respectively were 7.72 and 1.41 W. The front upper edge of the model is a quarter ellipse with a = 2.38W on the long half-axis and b = 1.39W on the short half-axis. The bottom edge has a circular transition with R = 0.113W. The gap between the bottom surface of the GTS model and the ground is h = 0.198W. The origin of the coordinates is located in the plane of symmetry between the bottom surface of the GTS and the back bottom surface of the model, shown in Figure 1(a). Four measurement points P1, P2, P3 and P4 were set on the base of the model, and the positions of the measurement points are shown in Figure 1(b). The distance of the measurement points in the span-wise direction is y = 0.8W, and the distance in the vertical direction is z = 0.8H. In this paper, eight various GTS models with Ra * = 2.0, 1.41, 1.0, 0.95, 0.87, 0.71, 0.56 and 0.5 are used to study the impact of Ra * on the asymmetric wake flow behind the GTS model. Here, the Ra * is defined as Ra * = H/W, where the Ra * value of the original model is 1.41. According to the previous studies (Grandemange et al., 2013;Zhang et al., 2022b), the Ra * is chosen to the range from 0.5-2, as listed in Table 1. It is worth noting  that GTS's height is changed to enable the variation of the aspect ratio. The computational domain with its corresponding boundary conditions used in this study are presented in Figure 2 (which has been presented in the study of Zhang et al., 2022b). Previous water channel test (McArthur et al., 2016) show that the wake structure behind the baseline GTS model is the equivalent between Re = 2.7×10 4 and Re = 2×10 6 , indicating a negligible influence of the Reynolds number on GTS's wake structure. Therefore, in order to save computational cost, all of the cases in this paper are simulated at Re = 2.7×10 4 . A uniform inlet velocity profile is used at the velocity-inlet, i.e. U inf = 6.16 m/s. To ensure that the boundary layer characteristics in this paper is consistent with that in the  Figure 2 is set as a rectangular no-slip region (stationary ground), the extra region has been set as moving wall. The detailed dimensions of the computation domain are depicted in Table 2 and Figure 2, which are detailed in the study of Zhang et al. (2022b). A commercial CFD software, Fluent 2020, was used in to investigate the turbulent flow based on IDDES with SST k-ω in the present study. The solver settings in this study are referred to the study of Zhang et al. (2022a) which investigated the effect of numerical methods on the wake structure. The momentum equation is used the bounded central differencing scheme (BCDS). The SIMPLEC is applied in the pressure-velocity coupling method. For the time derivative, the bounded secondorder implicit scheme is applied. The max iterations for each time step are set as 20, which give normalized values for all residuals lower than 10 −4 . The time step in this study is t = 0.01T inf . Here, T inf = W/U inf , defines a reference time for airflow traveling by one characteristic length (W) at incoming flow speed (U inf ). The data sampling of the flow field with a period of 241.6 T inf is initiated.

Mesh strategy and validation
The commercial grid generator Ansys ICEM-CFD is used to create a hexahedral mesh in the computational Figure 2. Diagram of the computational domain with its corresponding boundary conditions (Zhang et al., 2022b). domain. This mesh strategy has been widely adopted for the prediction of the turbulent flow around the GTS model (Rao et al., 2018a;Salari et al., 2004;Zhang et al., 2022a and2022b). The O-grid technique is applied to concentrate most of the elements close to the GTS model and in the wake region, as shown in Figure 3. The coarse, medium and fine meshes are used to investigate the mesh independence and the spatial resolution of the computational grid on the wake topology of the baseline GTS model. Figure 3(c)-(e) presents the grid distributions on the rear base of the baseline GTS model in the coarse, medium and fine grids. Note that the mesh distributions of three grids in the wall-normal direction (n + ) keep the same, but with the changes in the span-wise direction ( l + ) and stream-wise direction ( s + ). Here n + = u τ n/v, s + = u τ s/v and l + = u τ l/v, where u τ is the friction velocity, n is the distance between the first node and the GTS surface in the wall-normal direction, l and s are the GTS surface cell sizes in the stream-wise and span-wise directions, and v is the kinetic viscosity. Table 3 describes the main features of each utilized mesh. The n + mean is the mean value of the n + , and s + max and l + max are the maximum values of the l + and s + all over the surface of the GTS model. It should be noted that the O block is slightly adjusted in various aspect ratio cases to ensure the grid quality around the GTS model. Figure 4 compares the recirculation bubbles and the turbulent flow structures behind the GTS model predicted by the coarse mesh (first row), medium mesh (second row), and fine mesh (last row). The comparison of the time-averaged velocity streamlines is shown in Figure 4(a), (b) and (c). The mesh resolution has a great effect on the recirculation bubbles in the middle plane (y/W = 0). When the grid is too coarse, the wake flow topology is asymmetry, similar to Figure 3 (a) in the previous study of Rao et al. (2018a), but its asymmetry is weaker than that of Rao et al. (2018a), indicating the coarse IDDES cannot predict the shape of the recirculation bubbles. In contrast, the medium IDDES captures the main flow structures in the wake region as that predicted by the fine IDDES, with a sizable triangular vortex near the lower rear base of the GTS model and an elliptical vortex far away from the upper rear base (flow state I), which shows good agreement with that in the previous experiment of McArthur et al. (2016). This phenomenon also can be qualitatively explained by the comparison of the instantaneous wake structures. Moreover, Figure 4    grid refinement, i.e. the medium and fine meshes, more flow structures are resolved and captured, particularly in the separated flow regions. For a quantitative analysis of the influence of the mesh resolution on the wake flow of the GTS model, the averaged stream-wise (x-direction) and vertical (zdirection) velocity component profiles at eight different locations in the symmetrical plane (y/W = 0) are compared in Figure 5. The selected vertical lines are located at X/H = 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75 and 2.0. For the convenience of comparison, the X coordinate for the velocity profiles is obtained by the following relation: X = x+((u x,z /U inf )−1)×0.01. Here, x is the stream-wise coordinate, u x and u z are the values of velocity components in the stream-wise and vertical directions, respectively. Thus, the offset distance from its original vertical line represents the magnitude of velocity. The grey block in this figure indicates the GTS model, while the blue, black and red lines correspond to the coarse, medium and fine mesh, respectively. The general finding in Figure 5(a) is that the stream-wise velocity distribution in the wake region predicted by the medium IDDES shows good agreement with that of the fine IDDES, while the velocity profiles of the coarse IDDES are far from the medium and fine IDDES results, which is also confirmed by the vertical velocity component (u z ) profiles shown in Figure 5(b). This is because the coarse IDDES cannot accurately predict the shape and position of the recirculation bubbles behind the GTS model compared to those predicted by the medium and fine IDDESs. Furthermore, the flow topology of the coarse mesh is approximately anti-symmetric to those of the medium and fine meshes. From the above analysis of Figure 4, the flow state predicted by the coarse IDDES is similar to flow state II, while the flow topology results of the medium IDDES and fine IDDES are the same as flow state I.   Figure 6(a). The other flow state, which shows great agreement with the wake flow structure observed by the experiment of Schmidt et al. (2018) and the LES studies of Ortega et al. (2004) and Rao et al. (2018b), is named the flow state II, as visualized in Figure 6(b). Therefore, the validation of the dominant influence of the mesh resolution on the wake flow topology in the present IDDES simulations has been divided into two parts. The vortex cores predicted by the coarse IDDES are compared with the results of flow state I in Figure 6(c), while the medium and fine IDDES results are compared with the data of flow state II in Figure 6(d). Here, vortex A corresponds to the tiny vortex adjoining the large triangular-shaped vortex B. The vortex C indicates the smaller elliptical-shaped vortex. The saddle point in the vertical is denoted by S. Figure 6(c) compares the vortex cores predicted by the coarse IDDES (blue) that is similar to the flow state I, and the previous PANS and LES results obtained by Rao et al. (2018a) (red) and Rao et al. (2018b) (black). A significant finding is that the coarse IDDES misses the positions of the wake vortex cores and shows a large difference from the previous PANS and LES results, indicating an inadequate mesh resolution of the coarse mesh in this study. Figure 6(d) shows the comparison of the wake vortex cores results of the medium IDDES (blue) and fine IDDES (flow state I) (red) with the previous water channel test data measured by McArthur et al. (2016) (black). The general observation is that the medium IDDES and fine IDDES accurately capture the flow state II, and the wake vortex cores of the medium IDDES and fine IDDES show good agreement with the experimental data, especially for the vortex core positions of vortexes A and B. In  (Rao et al., 2018b) 10.78 × 10 6 0.5608 addition, the location of the saddle point S predicted by the medium IDDES is slightly closer to the experimental data than that of the fine IDDES case. Checking the core of vortex C in the medium IDDES, it's located a little far from the experimental data as compared to the fine IDDES. However, its distance between the origin point O and vortex core C only has a difference of 3.08%, when compared to that in the experiment. All indicate that the medium IDDES has a good prediction in capturing the turbulent wake topology of the GTS model. Last, a grid independence study is conducted to confirm the predicting accuracy of the GTS model's aerodynamic drag using IDDES. Table 4 compares the time-averaged aerodynamic drag coefficients predicted by the present IDDES and the previous LES. Here, C d is defined as C d = F d / (0.5ρU 2 inf S), in which F d is the timeaveraged drag force, ρ is the density of the air, and S is the area of GTS's rear base.
Considering the high accuracy of the resolved LES simulation in predicting the aerodynamic loads of the bluff body (Eric et al., 2013;Ortega et al., 2004;Salim & Ong, 2013), in Table 4, the C d value of the resolved LES is taken as the baseline. The drag coefficient predicted by the medium mesh with 5.10×10 6 elements is in good agreement with the fine mesh with 8.45×10 6 elements, while that predicted by the coarse mesh with 1.00×10 6 elements shows a larger discrepancy with the fine mesh. The relative error between the medium IDDES and the resolved LES (Rao et al., 2018b) is only 0.37%, proving a high accuracy of the medium IDDES simulation in the prediction of the aerodynamic drag of the GTS model.
All of results reported above indicate that the IDDES method works well for the prediction of the aerodynamic drag and the turbulent wake structures of the baseline GTS model. Furthermore, the overall good agreement with the water channel experiments (flow topology and vortex core position) and the resolved LES calculation (aerodynamic drag) allows us to select the present medium grid with the IDDES method to proceed in a more in-depth analysis of the effect of the aspect ratio on the wake flow topology behind the GTS model. Considering the high costs of the resolved LES simulation in predicting the wake flow structure (Rao et al., 2018b;Su et al., 2020), the medium IDDES is chosen for the following numerical simulations in this study.  Figure 7 shows the variation of the normalized timeaveraged drag coefficient (C d ) and drag coefficient area (C d S) with various Ra * . The length of the recirculation region both in the vertical plane and the horizontal plane is also presented in this Fig. Moreover, the normalized length of the recirculation region (L rh * in the horizontal plane and L rv * in the vertical plane) is normalized by the height of the GTS model (H * ). Both C d and C d S are normalized by the baseline C d and the

Flow topology calcification
The finding in Figure 7 is that the wake flow can be divided into three typical flow topologies (flow topology I, flow topology II, and flow topology III), according to the changing trend of C d value and the length relation between L rh * and L rv * . When the Ra * value is lower than 1.0, the length relationship of the recirculation region is L rh * > L rv * , and this flow topology is thereby named flow topology I (FT-I). As the Ra * value reaches 1.0, the characteristic length of the recirculation region in the horizontal middle plane is basically equal to that in the vertical middle plane (L rh * = L rv * ), which is defined as flow topology II (FT-II). For the flow topology III (FT-III) with the value of Ra * higher than 1.0, the characteristic length of the recirculation region in the horizontal middle plane is shorter than that in the vertical middle plane (L rh * < L rv * ). The changing trend of the GTS's C d value is significantly affected by its wake flow topology. When the wake flow switches into FT-I (Ra * < 1.0), the GTS's C d value decreases drastically with the increasing Ra * . As the Ra * value increases to Ra * = 1.0 (FT-II), the drag coefficient reaches the minimum value. As the flow topology falls into FT-III, the GTS's C d value increases obviously with the increasing Ra * . Furthermore, another finding in Figure 7 is that the normalized time-averaged aerodynamic drag coefficient area (C d S) presents a linear relationship with the value of Ra * , and the corresponding fitting equation is y = 0.7087x. The flow topology has a significant impact on the time-averaged drag coefficient while having a negligible influence on the drag coefficient area. It should be noted that, for the special flow topology (FT * , Ra * = 0.71), the interesting bi-stable wake flow phenomenon and the corresponding switching process between the wake flow states are observed in one numerical simulation, which will be detained analyzed in Section 3.2. Thus, only the results of the same flow state of FT * as that in FT-I are presented in Figures  8-15 in Section 3.1, in order to facilitate a more intuitive comparison.

Pressure distribution
In Figure 8, the visualization of the rear base pressure contours of the GTS models with various Ra * is presented. Here, the time-averaged pressure coefficient (C p ) is given by the equation: C p = (P-P ∞ )/(0.5ρU inf 2 ), in which P is the pressure on the rear base, P ∞ means the reference pressure, referring to the standard atmospheric pressure with the gauge pressure of 0 Pa in this paper (Deng et al., 2020). In this figure, the pressure contour shows a visible variation in various wake flow topologies. For the FT-III (Ra * = 2.0 and 1.41), the C p distribution near the bottom and that near the upper part has a large difference, while the C p distribution on the left side and right side of the rear base are basically symmetric. As the wake topology falls into FT-II, the C p distribution of the rear base is exactly symmetric in both span-wise and vertical directions. When wake topology falls into FT-I, the C p distribution becomes symmetric in the vertical direction, while being asymmetric in the span-wise direction.
For more quantitative comparisons, Figures 9 and 10 show the time-averaged C p profiles for the five vertical and horizontal sampling lines along the rear base, respectively. The positions of the sampling lines have been detailed in Figures 9 and 10(a). The vertical coordinates in all cases are normalized by H * (the GTS's height with various Ra * ) in Figures 9 and 10(b). The characteristic distributions of the C p profiles are determined by the flow topology to which it belongs. It can be seen that the C p values have a large difference in the higher and the lower sides along the vertical sampling lines for FT-III. The C p profiles in the FT-I and FT-II are both presented with less variation from the upper part to the lower part. The C p profiles along the span-wise direction in FT-1 shows large difference from left to right rear base, and the C p values on the left side of the rear base are significantly lower than the right side. However, this changing trend of C p profiles along the span-wise direction is disappeared immediately as the flow topology falls into FT-II. Figure 11 compares the iso-surface of the timeaveraged pressure coefficient C p = −0.16 around the GTS model in the front and top views. The C p isosurfaces are colored by the stream-wise velocity. As shown in Figure 11(a) and (b), the slant C p iso-surfaces of FT-III are found in the front view, while the C p isosurfaces in the top view is approximately paralleling the rear base, which is the reason why the C p distributions on the top of the vertical lines shown in Figure 9 have a large difference with the lower part. As the Ra * of the model decrease to the value of Ra * = 1.0 (FT-II), the isosurface is parallel to the GTS model both in the top view and front view. Furthermore, for the aspect ratio lower   Figures 8-10. Moreover, the various tilting directions of the iso-surface facilities to various flow structures, which will be further analyzed in Section 3.1.3.
The time histories of the gradient of the pressure coefficient in two directions and the results of the probability distribution function (PDF) have been shown in Figure 12 (the base C p gradient of FT * is presented in Section 3.2.1). The vertical C p gradient (∂C p /∂z * ) and the horizontal C p gradient (∂C p /∂y * ) are calculated using Eq. (1) and Eq. (2).

∂C
Where, the C p1 , C p2 , C p3 and C p4 are the C p signals collected using pressure probes P1, P2, P3 and P4 as shown in Figure 1(b). The results are both observed to exhibit the bimodal peak and the peak value has a large difference between each other, which confirms the various asymmetric characteristics of turbulent wake induced by the variation of the Ra * value (Grandemange et al., 2013). For FT-III, the vertical and horizontal PDF's peaks are observed at ∂C p /∂z * =0 and ∂C p /∂y * = 0, contributing to the asymmetric wake structures in the vertical plane and symmetric wake topology in the horizontal plane. In particular, the PDF's peak in FT-II occurs at the position of the C p gradient is 0, indicating that the wake structures are balanced along with the vertical and spanwise directions. As the Ra * decreases to Ra * < 1.0, the occurring positions of the vertical and horizontal PDF's peaks in FT-I occur with the value of ∂C p /∂z * = 0 and ∂C p /∂y * > 0, respectively, which means that the wake structures in FT-1 is symmetric and asymmetric along the vertical and span-wise directions. It can be concluded that in Figure 12 the wake topology of FT-III is asymmetric in the vertical plane but being symmetric in the horizontal plane, while the wake topology of FT-II is symmetric in both vertical and span-wise planes.
As Ra * = 1.0 decreases to Ra * < 1.0, the balanced wake   structures (FT-II) transform to be horizontal asymmetric and vertical symmetric characteristics (FT-I). Figure 13 compares the streamlines on the plane of y/W = 0 colored by the time-averaged stream-wise velocity. The value of the Ra * dramatically affects distribution of the stream-wise velocity in the near wake region, while the asymmetric characteristics on the plane of y/W = 0 are decided by the flow topology classification. In particular, the time-averaged streamlines are basically symmetric and the time-average velocity contours show symmetric in the vertical mid-plane (which is shown in Figure 13) when the wake flow topology belongs to FT-I and FT-II. As the wake flow structure falls into flow topology III (FT-III), a sizeable triangular vortex is presented on one side and a tiny elliptical vortex far is clearly visible on another side in the wake region, and this vertical asymmetric characteristic has a great agreement with the finding in the previous experiment (McArthur et al., 2016) and the LES simulation (Rao et al., 2018b). The aspect ratio (Ra * ) influence on the wake's asymmetric characteristics in the horizontal plane is shown in Figure 14. A pair of elliptical vortexes in the horizontal plane in FT-III and FT-II are presented. As the flow topology changes to FT-I, the visualization of the time-averaged streamlines on the plane of z/H * = 0 are asymmetric, it shows that a large triangular vortex is located near the left part and an elliptical tiny vortex is staying away from the right part of the wake region.

Asymmetric flow structure
The time-averaged flow structures for the cases with various Ra * in presented in Figure 15. The flow structures are represented by the mean iso-surfaces of Q = 3.1×10 3 s −2 . The length of the flow structures separated from the GTS's trailing edges is also highlighted. The flow structures are colored by the ω x * , which is defined as ω x * = ω x ·W/U inf , where ω x means the stream-wise vorticity. The general finding in Figure 11 is that the Ra * has a great influence on the wake flow topology. Three typical flow topologies can be recognized according to the size relationship of the wake flow separated from the trailing edges. For FT-I with the Ra * lower than 1.0, the length of the wake flow separated from the top trailing and the bottom trailing edge are basically the identical, while the right separated structure is obviously longer than that of the left separated wake structure. When the Ra * reaches Ra * = 1.0 (FT-II), the length of the lateral separated wake structures becomes the same, which is identical to the length of the upper and lower separated structures, and the wake topology falls into a balanced state in both of vertical and spanwise directions. As the Ra * value increases to higher than 1.0 (FT-III), the length of the lateral separated structures remains the same, while the length difference between the upper and lower separated structures becomes more significant, which leads to the symmetric recirculation bubbles on the plane of z/H * = 0.5 and asymmetrical recirculation bubbles on the plane of y/W = 0. Additionally, the configurations of the recirculation bubbles in the cases of Ra * = 1.41 and Ra * = 2.0 are quite different in the vertical mid-plane because the corresponding GTS's wakes stay in two different flow states of FT-III (Zhang et al., 2022b). Table 5 summarizes the Ra * effects on the near wake topology and asymmetric characteristics of the turbulent wake. In particular, when the Ra * is lower than 1.0, that the flow topology belongs to FT-I, and the wake structures are exactly symmetric on the plane of y/W = 0 and asymmetric on the plane of z/H * = 0.5. As the Ra * reaches to 1.0, that the wake flow topology falls into FT-II, and the wake structures get to a balanced state, showing the symmetric characteristics along with both vertical and span-wise directions. For the case of FT-III, the horizontal symmetric characteristics remain the same as that in FT-II, while the vertical asymmetric characteristics  appear immediately. Furthermore, three flow topologies significantly affect the changing trend of the drag coefficient and length relationship between horizontal and vertical recirculation region length. To sum up, the aspect ratio of the GTS model is the critical factor in the asymmetric characteristics in the near wake region, and the value of the Ra * = 1.0 is the inflection point of the GTS's asymmetry wake flow.

A random switching of the bi-stable flow states in the horizontal plane
The random switching of the bi-stable wake flow of the Ahmed body has successfully been observed in the previous numerical studies (Dalla et al., 2019;He et al., 2021). While for the GTS model (with a larger Ra * ), the bi-stable wake flow switching process has never been reported in the previous numerical simulation without any external condition drive, especially in the horizontal planes. Fortunately, during the numerical investigation of the Ra * effects on the topology of the GTS's wake flow in Section 3.1, a random switching of the bi-stable flow state in the horizontal plane in the case of Ra * = 0.71 (close to the crucial aspect ratio value of horizontal asymmetry of the Ahmed body with a square-back) is clearly observed. Thus, this section aims to investigate the random switching of the bi-stable asymmetric wake flow state in the horizontal plane. Owing to the dominant influence of the ground motion on the separation and evolution characteristics of the turbulent wake (He et al., 2022), the ground motion's effects on the corresponding random switching process between two adjacent wake flow state are also discussed in this section. It should be noted that the numerical simulation times for the moving ground (MG) and stationary ground (SG) cases are both increased to 497.2T inf from 362.4T inf , to accurately capture the random switching process. For the case of the MG condition, the no-slip region highlighted in Figure 2 is set to move at the same speed as U inf . Figure 16 shows the time histories ∂C p /∂y * with the corresponding probability distribution function results under the stationary ground (SG) condition and moving ground (MG) condition. For the SG case, two peaks of the PDF results located at two sides of the value of ∂C p /∂y * = 0 are visible, indicating two asymmetric mirrored flow states switching between each other appear in the horizontal plane in this numerical simulation. Furthermore, the higher value of Peak A than Peak B means that the time staying in flow state A (S A ) is longer than that in flow state B (S B ). This is because the special flow topology (FT * ) switches into S B from S A1 and falls into S A2 again. As shown in Figure 16(b), only one switching process from S A to S B is observed for the MG case during the same simulation period as the SG case. The result indicates that the ground motion has a great influence on the random switch between different flow state, and the ground motion is found to obviously inhibit the switching process and decrease the switching frequency of the bi-stable wake flow. Besides, Figure 16(b) presents that the ground motion contributes to higher Peak B than that of Peak A, suggesting a longer time of bi-stable wake flow staying in S B , indicating a significant effect of the ground motion on the distribution characteristics dominated flow state in the bi-stable wake flow. In particular, the bi-stable wake flow is dominated by S A in the SG case, while the corresponding dominated flow state becomes S B as the ground begins to move. Figure 17 shows the visualization of the aerodynamic pressure contours on the GTS's rear base with the flow states in the SG (SG-S A1 , SG-S B and SG-S A2 ) and MG (MG-S A and MG-S B ) cases. In Figure 17, the flow state has a great effect on the C p response at GTS's rear base. When the wake flow stays in S A1 and S A2 in the SG case, the C p value on the left side of the rear base is lower than the on the right side of the rear base, and the C p value near the upper and lower parts of the rear base is basically symmetric in the vertical middle plane, which is identical to the pressure distribution of FT-I shown in Figure 8. As the wake flow state is switched into flow state B in the SG case, the C p value on the rear base has the opposing characteristics in the horizontal plane as that in the S A1 and S A2 flow states. For the MG case, the C p distribution characteristics on the rear base in each flow state are similar to that in the SG case. The only difference is that the global C p value on the rear base in the MG case is relatively higher than that in the SG case. Figure 18 quantitatively confirms the comparison of the time-averaged C p profiles of the five lines in horizontal direction (the locations are detailed in Figure 18(a)) on the GTS's rear base with various ground motion conditions. Figure 18(b) shows that the distribution characteristics of the C p profiles of five flow states can be divided into two categories: (i) lower C p on the left base and higher C p on the right base (SG-S A1 , SG-S A2 and MG-S A ) and (ii) higher C p on the left base and lower C p on the right base (SG-S B and MG-S B ). Moreover, for the comparison of C p values on the right rear base between SG-S A1 , SG-S A2 and MG-S A , the MG-S A has higher C p values than that in the S A1 and S A2 flow state, owing to more replenishing energy near the right rear base which is caused by the higher speed reverse flow in the MG case presented in Figure 18. Similarly, the C p values on the left rear base in the MG-S B are obviously higher than that in the SG-S B , and this can be attributed to the higher speed reverse flow in the MG case shown in Figure 21. Compared to the SG case, the ground motion causes more airflow to flow reversely and impact the rear base in both two asymmetric mirrored flow states, which contributes to higher global C p value on the rear base.

Asymmetric flow structures
In this section, the turbulent flow structures in different flow states of FT * flow topology in the SG and MG conditions have been analyzed in detail. The time-averaged flow structures with the time-averaged iso-surfaces of Qcriterion at Q = 3.1×10 3 s −2 is presented in Figure 19. The length of the wake flow structures separated from the trailing edges in each flow state of the SG and MG cases is marked in this figure. The flow structures show that flow states significantly affect the flow structures in the wake region. When the wake flow is maintained in the flow state A (SG-S A1 , SG-S A2 and MG-S A ), the length of the wake flow structure separated from the left trailing edge is longer than the right side. For flow state B (SG-S B and MG-S B ), a longer flow structure is visible behind the left trailing edge. Furthermore, in both SG and MG cases, the switching process also changes the wake structures separated from the upper and lower trailing edges. For all flow states in the SG case, the wake structures separated from the upper trailing edge are always shorter than the lower side, while this relation becomes opposite in the MG case when the flow state switches into MG-S B from MG-S A . Additionally, the ground motion significantly increases the length of the wake structures separated from all trailing edges in the corresponding flow states compared to that in the SG case. Figure 20 shows the time-averaged streamlines in the horizontal mid-plane colored by the time-averaged stream-wise velocity in all flow states of the SG and MG cases. An interesting finding in Figure 20 is that the ground motion has a significant impact on the recirculation bubble configuration in S A and S B flow states of the FT * flow topology. A large triangular vortex near the left side while a tiny vortex keeping away from the right side are visible in the flow state A (SG-S A1 , SG-S A2 and MG-S A ) in both the SG and MG cases. While for the flow state B (SG-S B and MG-S B ) in the SG and MG cases, the asymmetric wake structures as the flow state A are clearly visible, showing the triangular vortex near the right side while a tiny elliptical vortex on the left side of the rear base. Furthermore, wake structures in S A are basically the same in both SG and MG cases. The only difference is the distribution of the reverse flow speed in the near wake region, and the MG case has a higher reverse flow speed near the right rear base and a lower reverse flow speed near the central rear base. Additionally, the asymmetric characteristics of the wake structures in MG-S B are more obvious than that in the SG-S B , owing to the longer remaining time of the flow state B in the MG case shown in Figure 16. The asymmetry characteristics of the wake structures between the MG-S B and SG-S B well explain the larger length difference between two lateral separated flow structures in the MG-S B than the SG-S B , as presented in Figure 19.
The time-averaged turbulence kinetic energy (TKE) contours in the horizontal mid-plane (z/H * = 0.5) in all flow states in the SG and MG cases are presented in Figure 21. Here, the TKE is normalized by U inf 2 . The distribution characteristics of TKE value vary obviously in the different flow states. In particular, for both SG and MG cases, the high-level TKE appears near the tiny elliptical vortex (right part of the near wake region) in the flow state A (SG-S A1 , SG-S A2 and MG-S A ). As the flow state switches into S B from S A in both SG and MG cases, the elliptical vortex occurs near the left edge of the GTS model shown in Figure 16, thereby resulting in a highlevel TKE distribution near the left edge. Furthermore, the ground motion has a significant influence on the TKE value in the flow state B while a negligible impact on the TKE distribution in the flow state A. Compared with the SG case, MG case shows a larger difference of the TKE value behind the left and right rear base, which confirms a more obvious asymmetric characteristics in MG-S B than SG-S B found in Figure 20.
For a more detailed analysis of the ground motion effects on the wake structures in the flow states A and B, Figure 22 and Table 6 compare the time-averaged stream-wise flow rate (R u ) on the rear cross-section. Four rear cross-sections are located at the rear base position and the normal height from the GTS's surface is the same as the ground clearance. The normalized time-averaged stream-wise flow rate is defined as R u * = R u /(U inf ·S * ), in which the S * is the area of the GTS's rear base with the case of Ra * = 0.71. The general observation is that for both the SG and MG cases, the R u * value of the left cross-section is higher than that of the right crosssection when the bi-stable flow state stays in S A (SG-S A1 , SG-S A2 and MG-S A ). As the bi-stable flow state in both SG and MG cases switches into S B , the R u * value of the right cross-section is higher than that of the left crosssection. The relation of the R u * value between the left and right cross-sections in S A and S B reasonably explains the wake structures in each flow state. This is because the increased flow rate causes a larger size of the shedding vortex separated from the corresponding GTS's trailing edge (Zhang et al., 2021). A stronger flow separation from  the trailing will lead to an increase in the rear base pressure and change the motion of the air in the near wake region, which causes the variations of the aerodynamic drag coefficients and dominant frequency in the near wake region (detailed in section 3.2.3). Furthermore, the switching process in both SG and MG cases also changes the R u * value of the upper and lower cross-sections, indicating a global impact of the wake flow state on the turbulent wake flow. Additionally, the ground motion significantly increases most of the stream-wise flow rate on the cross-sections, except for the R u * of the right crosssection in S A , which is the reason why the tiny elliptical vortex in the MG-S A is smaller than that in the SG-S A1 and SG-S A2 presented in Figure 20.  Tables 7 and 8 compare the time-averaged aerodynamic drag coefficients of all flow states in the SG and MG cases and the normalized time-averaged pressure coefficient integration (I Cp ) on the rear base of all flow states in the SG case, respectively. The I Cp is computed by summing the product of the area of each element by the value of the C p taken at the centroid of the element (nodal variables are interpolated at each cell centroid using cell shape blending or metric functions), for each element over the entire rear back. The final sum is then divided by the total area of the part, the equation is as follows.

Aerodynamic response
where, C pi is the time-averaged pressure coefficient taken at centroid of element i. A i is the area of element i. The time-averaged aerodynamic drag coefficients of SG-S A1 , SG-S B and SG-S A2 are 0.5364, 0.5417, and 0.5370, respectively. The C d value of SG-S A1 is in good agreement with that of SG-S A2 . There is a relative error  of C d between the SG-S A1 and SG-S B cases, revealing that the flow state greatly impacts the drag coefficients of the GTS model. When the flow state falls into SG-S B , the corresponding C d is larger than those of SG-S A1 and SG-S A2 cases. This is because the lower I Cp on the rear base of SG-S B (shown in Table 8) contributes to an increase in the aerodynamic drag coefficient during the vehicle operation, as that has been reported by Cheng et al. (2019) and Meile et al. (2016). For the case with moving ground, the time-averaged aerodynamic drag coefficients of MG-S A and MG-S B are 0.5419 and 0.5489, respectively. Compared to the aerodynamic drag coefficients in the SG case, . It confirms that the boundary layer influences the wake flow structure and the time-averaged force coefficients, and the accelerated airflow around the GTS model will cause an increase in the aerodynamic drag coefficient.
To assess the dominant dynamics of the wake structure, the comparison of the dominant frequency on the near wake region of all flow states in the FT * flow topology is performed. The Fast Fourier Transform (FFT) algorithm is employed to obtain the dominant frequency based on the time history of the span-wise velocity component. The probe for the collection of the span-wise velocity of flow state A (SG-S A1 , SG-S A2 and MG-S A ) and flow state B (SG-S B and MG-S B ) is located at the position of X/H * = 1.5 (stream-wise direction), Y/H * = 0 (span-wise direction) and Z/H * = 0.5 (vertical direction). Here, the Strouhal number (St) is definded as St = f ·H/U inf , f is the frequency, and the power spectral density (PSD) in this figure is normalized by the corresponding maximum value. Figure 23 shows that the dominant St of the asymmetric mirrored flow states is quite different in both SG and MG conditions. In particular, the dominant frequency in SG-S A1 and SG-S A2 are both St = 0.104, while it increases to St = 0.166 when the flow state switches into SG-S B from SG-S A1 . This reveals that the dominant frequency is highly affected by the belonging flow states; even the adjacent flow states are observed for the same ground condition. Furthermore, the ground motion obviously impacts the dominant shedding frequencies of both flow states A and B. Compared to the SG case, the dominant shedding frequencies of the asymmetric mirrored flow states A and B increase to St = 0.187 and St = 0.145 in the MG case, respectively. This might be because the ground motion drastically increases the underbody flow rate (presented in Figure 22 and Table  6), and the complement of the underbody flow thereby increases the corresponding response frequency in the wake region (e.g. SG-S A and MG-S A , SG-S B and MG-S B ).

Conclusion
In this paper, the aerodynamic performance, wake topology and asymmetric characteristics of the GTS models with Ra * = 0.5, 0.56, 0.71, 0.87, 0.95, 1.0, 1.41 and 2.0 have been investigated. The Ra * effects on the GTS's wake flow structure and asymmetric characteristics in the near wake region have been studied using IDDES simulation. The influence of the ground motion on the flow characteristics in the asymmetric mirrored bi-stable flow states A and B in the special flow topology Ra * = 0.71 has been analyzed. The primary conclusions are summarized as follows.
(1) Three typical flow topologies (FT-I, FT-II and FT-III) are classified according to the changing trend of the drag coefficient and the length relation between the vertical and horizontal recirculation region. In the FT-I (Ra * < 1.0), the horizontal recirculation region length is longer than the vertical, and the drag coefficient decreases drastically with the increase of Ra * . When the wake flow stays at FT-II (Ra * = 1.0), the wake flow structures get into a balanced state and the drag coefficient reaches the minimum value. When wake flow keeps into FT-III (Ra * > 1.0), the relationship between the horizontal and vertical recirculation region length and the changing trend of the drag coefficient with the aspect ratio value are opposite to the FT-I. (2) The vertical and horizontal asymmetric characteristics of the GTS's wake structures can be transformed between two adjacent flow topologies with the variation of the aspect ratio. When the flow topology belongs to FT-III, the wake flow structures are asymmetric in the vertical mid-plane and symmetric in the horizontal mid-plane. As the flow topology falls into FT-II, the wake flow remains symmetric along with the span-wise direction, and the wake structure becomes symmetric in the vertical mid-plane. As the flow topology falls into FI-I, the vertical symmetric characteristics remain the same as that in FT-II, and the horizontal asymmetric characteristics appear simultaneously. (3) For the case of Ra * = 0.71 of the GTS model, the wake asymmetric characteristics are exceedingly active, and a random switching process between two asymmetric mirrored flow states S A and S B in the horizontal plane is successfully observed. The dominant shedding frequencies of the asymmetric mirrored flow states are quite different from each other. The ground motion does not affect the bistable flow states category while reversing the dominated flow state and inhibiting the corresponding switching process. Furthermore, the ground motion drastically enhances the asymmetric characteristics of the wake structures in the S B . In addition, the ground motion also significantly affects the pressure distribution, flow configuration, TKE level, flow rate and dominant shedding frequency in the near wake region in each flow state.
The current research revealed mechanism of the aspect ratio on the GTS's asymmetrical wake flow characteristics and the corresponding aerodynamics performance. Due to the limitation of the simulation resource, a detailed wake evolution analysis is not presented in this paper. Therefore, future work will perform the long-time numerical simulations and the POD (proper othorogonal decomposition) analysis will be carried out to further investigate the motion modes of the GTS's asymmetrical wake flow and reveal the switching mechanism of bi-stable flows.