Numerical Analysis of the Flow around Two Square Cylinders in a Tandem Arrangement with Di ﬀ erent Spacing Ratios Based on POD and DMD Methods

: To more clearly understand the changes in ﬂow characteristics around two square cylinders with di ﬀ erent spacing ratios, the main mode of the ﬂow ﬁeld was extracted by using the Proper Orthogonal Decomposition (POD) and Dynamic Mode Decomposition (DMD) methods. The changes in the main mode of the ﬂow ﬁeld at di ﬀ erent spacing ratios and the di ﬀ erence of the time series were analyzed and compared. This processing can separate the mixed information in the ﬂow ﬁeld and obtain the dominant modes in the ﬂow ﬁeld. These main modes can clearly reﬂect the dominant ﬂow characteristics in the ﬂow ﬁeld. The analysis results show that when L / D = 2, the ﬂow ﬁeld structure is consistent with the ﬂow ﬁeld around a single square cylinder. When L / D = 2.5–3.5, the vortex shedding from upstream cylinders combines with the vortex near the downstream cylinders. This mutual coupling causes a signiﬁcant change in the drag coe ﬃ cient value of the downstream cylinder. When L / D = 4, the main vortex from the upstream cylinder can be completely shed, which means that the upstream and downstream square cylinder vortices start to become independent. The main focus of this paper is to use the advantages of POD and DMD to obtain several modes with higher energy in the ﬂow ﬁeld. Furthermore, it can be considered that these main modes can fully reﬂect the ﬂow characteristics of the ﬂow ﬁeld.


Introduction
Square cross-section structures commonly exist in the fields of civil and marine engineering. The problem of flow around square cylinders has been a popular topic in structural wind and ocean engineering research in recent years. In practice, structures rarely exist in isolation and more often present a certain form of arrangement. Two cylinders in a tandem arrangement are one of the simple and typical distribution forms, which, to some extent, represents the complexity of structural interference in the flow field.
Liu [1] studied the aerodynamic changes of the flow around two square cylinders with different spacing ratios under several Reynolds numbers (Re) by means of experimental methods. Sakamoto [2] carried out experiments under Re = 2.76 × 10 4 and Re = 5.52 × 10 4 , and obtained three spacing ratio intervals for the aerodynamic change of the upper and lower square cylinder. Kim [3] examined the test flow around two square cylinders with Re = 5.3 × 10 3 and Re = 1.6 × 10 4 , and found that the aerodynamic coefficient curve changed significantly at the point where the spacing ratio was 2.0 and 2.5, respectively. characteristics of the flow field is only by means of experiments or CFD. This means the data obtained from these methods generally contain a variety of information components, such as various scales of vortex display, on one image, and it is difficult to identify the contribution of different vortex structures to the flow field based on these instantaneous flow characteristics. The main focus of this paper is to use the advantages of POD and DMD to obtain several modes with higher energy in the flow field. Furthermore, it can be considered that these main modes can fully reflect the flow characteristics of the flow field. Based on the characteristics of these main modes, comparison and analysis of the changes in the flow field around two cylinders in a tandem arrangement with different spacing ratios at Re = 1.6 × 10 4 is undertaken. It should be noted that the source data of the analysis are based on CFD simulations.

CFD Simulation
A quadrilateral grid was used for calculation. The grid division and area size are shown in Figure 1, where D is the side length of the square cylinder and L is the distance between the square cylinders' centroids. The mesh of the single cylinder is similar to that of two cylinders. The turbulence model chosen for calculation is the RNG k-ε model. The no-slip boundary condition is imposed at the cylinder's surface, the left semicircle boundary as the velocity inlet condition, and the right semicircle as the pressure outlet. Using the semi-implicit method for pressure-linked equations (SIMPLE) algorithm, the first order upwind is used as the format of temporal discretization, and the spatial discretization precision is second order. The non-dimensional time-step size is ∆s = ∆tU/D= 0.0002 × 3.6/0.065 = 0.0111, including more than 500 time-steps in each period of vortex shedding [30], and performing 20 inner iterations within each physical time-step. Simulations were performed using ANSYS Fluent. Corresponding grid independence test was carried out, taking the two square cylinders with a spacing ratio of 2.5 as an example, and the results are shown in Table 1. Finally, the medium mesh was selected for calculation. The size of the first layer of the grid near the wall is y 1 (∆y = y 1 /D= 0.0003); calculation results show that y + max ≈ 1. The lift coefficient RMS value C l ' and drag coefficient average value C d are shown in Figure 2. The simulation results of the upstream cylinder are very close to the test values of Liu [1]. However, as can be seen from the aerodynamic coefficient variation trend and transverse velocity U distribution (Figure 3), when L/D equals 2, the airflow circumvents the upstream square cylinder and then reattaches to the downstream cylinder, and a complete vortex shedding is then formed behind the downstream square cylinder. The wake flow pattern is very similar to the flow around the single square cylinder ( Figure 4). When the spacing ratio changed from 2 to 2.5, the negative C d value of the downstream square cylinder changed to positive. The vortex shedding from the upstream square cylinder failed to reattach to the downstream cylinder on both sides in time. Instead, it crashed into the windward side of the downstream square column, causing the flow around the downstream square column more complex. As shown in Figure 3, the wake flow pattern of L/D = 2 is quite different than that of the others; this spacing ratio fits with the turning point of the force coefficient in the test of Kim [3]. As the spacing ratio increases to 4, it can be seen that there is a relatively complete vortex sandwiched between the upper and lower cylinder, which is generated by the airflow passing around the upstream square cylinder. The direct results obtained by CFD simulation can only show the change of a variable value over time in a single reaction flow field, and these data contain both temporal and spatial information. In this work, POD and DMD methods are used to separate the main modes of the flow field and the corresponding time coefficient. As described above, with this processing method, not only can instantaneous images of the flow field be obtained, but flow field modes with a clear energy proportion (degree of contribution) can also be obtained. The calculation results of the flow around the single square cylinder are shown in Table 2, and the results are in good agreement with the test values. It can be seen from the comparison between CFD simulation results and experimental values that the flow field obtained by calculation has high reliability, and the flow field data obtained by calculation can be used for the following analysis.

Proper Orthogonal Decomposition
POD has a complete mathematical description in mathematics. Its essence is to find a set of optimal basis functions in the space domain to express a known function. It is essentially a linear decomposition [31]. After a series of mathematical descriptions and transformations, generally in actual use, POD can be expressed as shown in Equation (1), where ϕ m (x) is the POD mode and a k m is the corresponding coefficient of the mode. For rigorous mathematical derivation, please refer to Lumley [13].
In the snapshot POD method proposed by Sirovich [32], M represents the number of snapshot samples. The modes and features identified by POD are ranked based on energy. Therefore, the first few order modes obtained by the POD method reflect the main flow characteristics of the flow field, and usually the dominant character of flow fields can be reconstructed well using these modes. In fact, the POD method is basically consistent with the singular value decomposition (SVD) and principal component analysis (PCA) methods [25]. The SVD formula is shown in Equation (2), where W is the variable matrix of interest, and each column of vectors represents a time sample. If it is a two-dimensional spatial distribution, it should be rearranged into a one-dimensional vector. U returns the POD mode matrix, V returns the time series corresponding to each order mode, and S is the energy ratio of each order mode, where * represents the conjugate transpose.

Dynamic Mode Decomposition
We can redistribute the variable matrix W into W A and W B , where W A contains the first N-1 snapshot, and W B is the last N-1 snapshot (N is the total number of samples): Suppose W B can be obtained by W A transformation, as shown in Equation (4), and the singular value vector of W A is decomposed into W A = USV * . The solution of the conversion matrix F can be further simplified to Equation (5) [25], and then the eigenvalue solution of the matrix F is derived to obtain the eigenvalue µ j and the eigenvector Λ j .
The frequency corresponding to each mode is shown in Equation (6), where f s is the snapshot sampling frequency. The modal flow growth rate (decay rate) can be expressed as Equation (7).
Finally, the DMD mode can be expressed by Equation (8), where Λ j is the eigenvector corresponding to the mode.
The ordering of DMD modes also requires reflecting the degree of contribution of this order mode to the flow field but, at present, this sorting method has not been unified like POD. Here, several sorting methods exist: sorting according to the modal 2-norm [24,33]; sorting rules developed by Kou [34]; and the sorting method used in Hu's paper [35], as shown in Equation (9). The number of snapshot samples used in decomposition is 500, and the sampling frequency of snapshot samples is about 714 Hz. POD decomposition is performed on the transverse velocity U which is calculated by CFD.

Results and Discussion
The POD decomposition of the sample matrix W can be performed directly, so that the first-order POD mode obtained is the average flow field. It is also possible to subtract the mean value of the snapshot samples from each value in the matrix W, so that the obtained mode does not include the first-order mode of the average flow field. In order to compare with the DMD method, this article adopts the former method for POD processing.
The energy distribution of the main POD modes of the flow field around a square cylinder and the flow field around two square cylinders in a tandem arrangement at different spacing ratios is shown in Figure 5. The first-order mode occupies the dominant energy, which reflects the character of the average flow field or static flow field, as shown in Figure 6. It can be seen that modes 2 and 3, 4 and 5, 6 and 7 actually exist in pairs, from their corresponding time coefficients, it can also be seen that a phase difference is formed between the two modes, and the frequency distribution obtained after the spectral transformation is also consistent, which reflects the alternating shedding of the vortex, which is consistent with Wang's related research [36]. In this article, the operating conditions, with the exception of L/D = 2, do not list another mode of lag. At the same time, it can be seen from Figure 6 that each order POD mode of L/D = 2 is very similar to the flow around a single square cylinder (Figure 7). The order 2-3 is an anti-symmetric structure, and the order 4-5 is a symmetric structure. The vortex structure changes from large to small as the order of modes increases, corresponding to the gradual decrease in vortex energy. Except for the average flow field, the main frequency obtained by the second to third order mode spectrum transformation with the largest energy ratio is very close to the main frequency obtained by the upstream and downstream square cylinder lift coefficient spectral transformation, and dominates the pulsating field. Extracting the first few modes for reconstruction, it can be seen that the reconstructed images of the order 1-5 and the order 1-7 are very consistent with the original image, which shows that the order 1-5 POD modes can reproduce most of the characteristics of the flow field. As the modulus increases, the accumulated modal energy changes little, resulting in little change in the reconstructed image. As shown in Figure 8. When the spacing ratio between the upstream and downstream square cylinder changes from 2 to 2.5, the POD mode is no longer similar to the flow field around the single square cylinder. The wake area of each mode expands, and the vortex core becomes smaller. The character of the flow field is quite different from that of L/D = 2, but it can still be seen that the third and fifth mode show obvious asymmetrical and symmetrical structures, respectively. When L/D equals 3 or 4, the symmetry of the POD mode is weakened, which may be due to the fact that these POD modes are simultaneously doped with eddy forms of other frequencies. As shown in Figure 9, when L/D = 2.5, the spectral transformation corresponding to the time coefficient of the POD mode has only a single main frequency, but when L/D > 3, some modes have multiple main frequencies. Table 3 lists the main frequencies obtained by the spectral transformation of the main POD mode time coefficients with different spacing ratios, and the main frequency obtained by the time series spectral transformation of the upstream and downstream square cylinder lift coefficients. From these frequency distributions, it can be seen that when L/D equals 2 or 2.5, the dominant frequencies of the first and third modes are the eddy frequency or its multiplied frequency. When the spacing ratio is gradually increased, the second and third-order mode frequencies are basically unchanged, but the frequency of the fourth to seventh modes gradually increases. When the spacing ratio L/D = 4, the frequency of the fourth to seventh mode is close to the multiple relationships with the vortex shedding frequency of the square cylinder. From the mode diagram of L/D = 4, it can be seen that the main vortex of the upstream and downstream square cylinders begins to separate from each other, that is, after the airflow flows through the upstream square cylinder, the vortex shedding that dominates the energy can be formed relatively completely, and then flap on the downstream square cylinder.     As described above, the mode obtained by the POD method may contain multiple frequencies. In practice, in civil engineering or marine engineering structures, the vibration response of the structure is often closely related to a certain frequency, such as vortex-induced vibration (VIV). Therefore, if we can obtain the frequencies with relatively high energy in the flow field and the flow characteristics corresponding to these frequencies one-to-one, it will be beneficial to our analysis of the structural vibration mechanism. The dynamic modes obtained through DMD can display such information well. In this study, DMD processing was also performed on the transverse velocity U, and then sorted according to the energy occupied by each mode. These DMD modes appear in pairs, with opposite frequencies. The DMD mode energy at different spacing ratios is shown in Figure 10. In most papers, the first-order mode obtained after sorting by mode amplitude or mode energy is the average flow field, that is, the 0 Hz mode [33][34][35]. However, after sorting according to the selected sorting rules in this paper, the first-order DMD mode obtained under the condition of L/D equals 2.5, 3, or 3.5 is not the 0 Hz mode, but is close to or at the same order of magnitude as the 0 Hz mode in terms of energy proportion. Such analysis results also appear in the analysis of Yu [26].
The distribution of the Ritz values of the DMD modes at each spacing ratio is shown in Figure 11, where the small green circles represent on or inside the unit circle, and the red solid squares represents outside the unit circle. It can be seen that most of the modes at each spacing ratio are in a stable state or a periodic state. When L/D = 2.5, there are two frequencies with a large growth rate around 0 Hz, but their energy ratio is very small, so they can be ignored. Table 4 lists the growth rate (decay rate) of DMD modes with different spacing ratios, and the frequency and structure of the modes are shown in Figure 12. When L/D = 2, the first-order mode is 0 Hz and occupies the dominant energy. The attenuation rate of the first few modes is close to 0, which represents a stable state. Their frequency is basically the same as the POD mode, and close to the vortex frequency or a multiple of it. When L/D = 2.5, it can be seen that there are multiple modes with larger energy, each of which representing a different frequency, meaning that there are many main frequency components in the flow field. The structure of the flow field is more complicated, but the frequency of the first, second, fourth, and fifth orders decays quickly, and the corresponding frequency and the vortex frequency deviate slightly. These characteristics reflect the interaction between the front and rear square cylinders in the whole flow field well. When L/D = 3, the energy ratio of the static flow field is close to other sixth-order energies, and in a stable state. These modes together dominate the flow field. At this time, the coupling of the upstream and downstream vortices is relatively weakened. When L/D = 4, the average flow field returns to the first order, and the energy ratio of the first mode is very large. The energy distribution in the flow field is similar to L/D = 2, indicating that the number of main frequency components in the flow field is small, and the flow field structure is clear. The coupling degree of the vortex shedding from the upstream cylinder and the vortex around the downstream cylinders is weakened. It can be seen that the results obtained by DMD analysis and POD analysis are in agreement.

Conclusions
Through the CFD simulation, the data of the flow field around two cylinders in a tandem arrangement in different spacing ratios and the flow field around a single cylinder were obtained. These data were then compared with experimental data to verify the reliability of the simulation. Finally, POD and DMD analysis was performed on these data, and a comparison was made of the main modal features extracted to fully understand the change of flow field characteristics at different spacing ratios. The following conclusions can be drawn: (1) In the POD analysis, when the spacing ratio L/D = 2, the POD mode is similar to the flow around the single square column, and the frequency of the first few modes is close to the frequency of the vortex shedding from the square cylinder or their multiples. When the spacing ratio increases to 2.5, POD modes are greatly different from L/D = 2, and the front and rear square column vortexes are coupled with each other, with an obvious interference effect. When the spacing ratio is 4, it can be seen from the POD mode whose energy proportion is only second to the average flow field that the vortex shedding that dominates the energy can be formed relatively completely after the airflow flows through the upstream square cylinder.
(2) In DMD analysis, when the spacing ratio L/D = 2, the DMD mode and frequency are consistent with the POD mode, and the main modes are in a stable state; when L/D equals 2.5 or 3, there are multiple main modes with relatively high energy, indicating that the flow field is composed of multiple main frequency components; when L/D = 4, the DMD mode energy ratio is similar to L/D = 2, and the main mode frequency with a large growth rate is close to the frequency of vortex shedding or its multiples, indicating that, at this spacing ratio, the upstream and downstream square cylinder vortices start to become independent.
(3) According to the above analysis results, it can also be seen that POD and DMD achieve the separation of complex information in the flow field. It can be seen from the text that the lateral velocity image calculated by CFD is difficult to clearly describe the flow field characteristics caused by the spacing ratio. A series of modes sorted by energy, obtained by POD decomposition, represent eddy structures with different participation degrees mixed in the flow field. By comparing the first few modes with the higher energy, we can clearly see the change of the flow field eddy caused by the change of the spacing ratio. When the flow field itself is more complex, the POD mode may be mixed with multiple frequency components. The mode obtained by using DMD directly corresponds to a specific frequency, so that we can obtain the dominant frequency of the flow field and the corresponding flow structure more clearly. The increase in frequency has brought benefits for us to further study the flow-induced vibration of the structure. In future work, assuming that the accuracy of CFD simulations of complex cross-sections such as bridge cross-sections can be improved, and that such simulations can be performed more cheaply than at present, the complex flow fields around these sections can be decomposed to further understand the flow patterns, and even the wind-induced vibration mechanism, in these cases.

Acknowledgments:
The authors thank the editors and reviewers whose constructive suggestions led to an improved paper, which are gratefully acknowledged. height of the first layer near the wall (m) y +max max height of the dimensionless first layer near the wall µ j eigenvalue of DMD Λ j the eigenvector of DMD ω k (x)

Conflicts of
Snapshot data matrix ϕ m (x) POD modes Φ j DMD modes Φ j energy of each DMD modes