1. Introduction
In the domain of computational fluid dynamics (CFD), the pursuit of simulating complex fluid flows with high accuracy and efficiency has been a driving force behind technological advancements in various industries, from aerospace and automotive engineering to environmental science and biomedical research. However, fine discrete computational meshes are often used to develop accurate numerical solutions for CFD problems using high-performance computing. Accurate numerical solutions usually require higher-order discretization schemes for complex partial differential equations, which are prone to produce oscillations [Reference Hannoun and Alexiades10].
An important technique that has revolutionized the capabilities of complex CFD simulations, to improve the efficiency of flow solvers and the accuracy of numerical solutions, is adaptive mesh refinement (AMR). AMR is a computational approach that enables users to optimize their simulations by dynamically adjusting the resolution of the computational mesh in response to evolving flow features and phenomena, and reducing the computational cost [Reference Li15].
Early studies in the field of adaptive mesh refinement for CFD problems can be traced to the 1980s and 1990s, including the works of Berger and Oliger [Reference Berger and Oliger4], Bell et al. [Reference Bell, Berger, Saltzman and Welcome2], Friedel et al. [Reference Friedel, Grauer and Marliani9], and Berger and Leveque [Reference Berger and Leveque3]. The development of AMR has continued to evolve with advancements in numerical methods, computer hardware and software tools, and it has been applied to a variety of fluid dynamics and computational science problems [Reference Antepara, Balcázar and Oliva1].
Li [Reference Li13] developed a two-dimensional (2D) AMR method derived from the qualitative theory of differential equations. The AMR method refines a given mesh based on the numerically computed velocity fields. The efficiency and accuracy of the AMR method has been verified using the accurate locations of singular points, asymptotic lines and closed streamlines [Reference Li12], and against widely used CFD benchmark experiments including the lid-driven cavity flow [Reference Lal and Li11], the 2D unsteady flow past a square cylinder [Reference Li14], the backward-facing step flow [Reference Li and Li17] and 2D flow over a wall-mounted plate [Reference Li and Lal16]. In particular, the AMR method has been shown to be useful for capturing localized flow features such as accurate location of the centre of vortices within the refined cells [Reference Lal and Li11, Reference Li and Lal16]. The AMR method is robust [Reference Lal and Li11], low-cost [Reference Li15] and can be applied to any incompressible fluid flow [Reference Li, Mallinson, Yangbo, Kaoru, Cluckie and Smedt18]. The previous works, for example, [Reference Li14, Reference Li and Li17], considered the accuracy of the 2D AMR method with two refinements and used the finite volume methods. We showed that the twice-refined cells contain the centre of vortices. Since the 2D AMR can be applied a finite number of times, more accurate centre of vortices can be obtained if the numerical velocity fields are accurate enough and more refinements are carried out. The accuracy of the AMR method depends only on the accuracy of the numerical methods used.
Flow around a square cylinder is a classic CFD problem and has been studied by many researchers (see, for example, [Reference Dhiman, Chhabra and Eswaran7, Reference Erturk and Gokcol8, Reference Perumal, Kumar and Dass19]). The flow pattern characteristics are related to factors such as the Reynolds number ( $Re$ ) of fluid, the dimension of the square cylinder and the space size. The complex flow characteristics in the vicinity of a square cylinder have always been a focus area and a challenge in academic research. The outcomes of this study can be applied to problems in aerodynamics, wind engineering and electronics cooling. One of the examples is the design of a building. The flow around a building has high- and low-pressure areas. These areas significantly affect building structural integrity, energy efficiency and indoor comfort.
It is noted in the literature that the flow around a square cylinder is steady for $5\le Re < 60$ and a fixed blockage ratio of $1/8$ [Reference Dhiman, Chhabra and Eswaran7, Reference Perumal, Kumar and Dass19]. Accordingly, the present study aims to use the AMR method of Li [Reference Li13] and provide estimates of coordinates of the locations of the centre of vortices for steady incompressible flow past a square cylinder. The AMR method is tested on eight cases considering flows with different Reynolds numbers, that is, $5\le Re\le 50$ . Hence, the main contribution of our work is showing the accuracy of the AMR method using four refinements and the finite-element method. Furthermore, we show that the AMR method can capture the location of the centre of vortices within the fourth refined cell. The efficiency of the AMR method in identifying accurate locations of the centre of vortices is validated against the available numerical results found in the literature in which a much finer mesh was used for the study. The methods in the literature investigating flow around a square cylinder (for example, [Reference Erturk and Gokcol8]) are limited to steady incompressible flow only. However, the AMR used in this paper can be applied to both steady and unsteady incompressible flows.
The remainder of this paper is structured as follows. Section 2 outlines the methodology, presenting the governing equations, the computational mesh and the flow solver. The numerical results of the eight test cases and the accuracy verification of the 2D AMR method are presented and discussed in Section 3. The conclusion follows in Section 4.
2. Methodology
This section briefly describes the governing equations, the computational domain and the mesh structure, and the flow solver.
2.1. Governing equations
We consider the finite-element discretizations of the 2D unsteady, incompressible Navier–Stokes equations [Reference Cannon and Knightly5], which contain the continuity equation and the momentum equations in two directions defined by
where $\mathbf {V}=(u, v)$ denotes the velocity field in 2D with u and v as the velocity components in the x- and y-directions, respectively, $\nu $ is the kinematic viscosity, $\rho $ is the fluid density, and P represents the scalar pressure.
2.2. Computational domain and boundary conditions
The geometry of the computational domain was set as a 2D rectangular channel with length L and width H, as shown in Figure 1. The two ends of the channel were the inlet and outlet. A square cylinder with a side length of $D=1$ is placed in the middle of the channel centred on the y-axis such that the blockage ratio of the channel, $D/H$ , equals $1/8$ . To ensure that the distance from the computational inlet to the cylinder does not influence the accuracy of the numerical solution, an essential inlet location of approximately 10 has been reported [Reference Sohankar, Norberg and Davidson20]. Thus, we set the inlet length $L_1=10$ . The channel length L is set to 30.
The boundary conditions used are described in Figure 1. Following [Reference Erturk and Gokcol8], the inlet velocity field was specified as a parallel flow with a parabolic horizontal component given by $u=U(y)= 1-(y-4)^2/16$ . The outlet was set as the free outflow boundary ( $P=0$ ). The channel’s top and bottom boundaries and the square cylinder walls are imposed as rigid surfaces with the no-slip condition ( $u=v=0$ ). The Reynolds number is considered the main parameter that changes the flow behaviour and is defined as $Re=U_{\max }D/\nu $ , where $U_{\max }=1$ and $\nu $ is the kinematic coefficient of the viscosity [Reference Erturk and Gokcol8, Reference Perumal, Kumar and Dass19]. Eight test cases for flows with different Reynolds numbers ${(Re =5, 10, 15, 20, 25, 30, 40}$ and 50) are considered in this study.
2.3. Computational mesh and the flow solver
For all test cases involving different Reynolds numbers, the spatial step size in both the x and y directions is $1/4$ . The initial mesh ( $\text {Mesh}_0$ ) has 3824 quadrilateral cells (size $1/4\times 1/4$ ) with 3984 nodes. Figure 2 shows the initial mesh $\text {Mesh}_0$ . The Finite Element Analysis simulation Toolbox (FEATool, version 1.16.3, https://www.featool.com/) in MATLAB (R2023a) was used to numerically solve the unsteady 2D Navier–Stokes equations [Reference Cannon and Knightly5] for the eight test cases.
This study aims to verify the accuracy of the existing adaptive mesh refinement method with accurate numerical velocity fields. The time-dependent Navier–Stokes equations are solved using FEATool until both velocity fields and pressures converge to a steady state. The highest-order scheme implemented in FEATool is a second-order implicit Crank–Nicolson time-stepping scheme [Reference Crank and Nicolson6], which was used for the time-dependent settings. We assumed that the convergence is achieved when the numerical solutions (velocity fields) are computed with residuals smaller than $10^{-10}$ for both u and v. Thus, the time-dependent simulations terminated when either the central processing unit (CPU) time reached the specified simulation time or the normed changes in the dependent solution variables were smaller than the time-stopping criteria (tolerance = $10^{-10}$ ), and at the same time, we make sure velocity fields and pressure show the steady state. Therefore, in this work, only steady-state solutions are considered.
3. Numerical results and accuracy verification
In the present study, we apply the same AMR method from Li and Lal [Reference Li and Lal16]. We refer the reader to the work of Li and Lal [Reference Li and Lal16, Section 3] for a detailed description of the method. A cell-by-cell AMR approach is followed on the quadrilateral meshes. A cell is refined by connecting the mid-points of opposite sides into four equally smaller quadrilateral cells. Applying the AMR method to the initial mesh, $\text {Mesh}_0$ , produces the first refined mesh, $\text {Mesh}_1$ , and by repeating the procedure three more times, we obtain the fourth refined mesh, $\text {Mesh}_4$ . The coordinates of the centre of the isolated refined cells in the area of interest are then used to estimate the centre of vortices.
The results of the eight test cases are presented in this section. We verify the accuracy of the AMR method by comparing it with the reference centre of vortices, $(x_{\text {ref}}, y_{\text {ref}})$ , presented by Erturk and Gokcol [Reference Erturk and Gokcol8]. There, the authors used the streamfunction and vorticity formulation of the Navier–Stokes equations and employed a very fine mesh with close to 1 040 000 nodes.
3.1. Refined meshes
Figure 3 shows the fourth refined mesh ( $\text {Mesh}_4$ with 13 223 nodes) for the flow with $Re = 5$ . The mesh is refined mostly near the channel’s top and bottom boundaries and near the square cylinder.
Figures 4(a)–4(e) show the zoomed-in sections of $\text {Mesh}_0$ to $\text {Mesh}_4$ in the vicinity of the square cylinder for the flow with $Re = 5$ . In $\text {Mesh}_0$ , a cell to the right of the square cylinder is marked with a red square. The first refinement of the cell marked with a red square is shown in the once-refined mesh, $\text {Mesh}_1$ . Similarly, the cell’s second, third and fourth refinements are shown in $\text {Mesh}_2$ , $\text {Mesh}_3$ and $\text {Mesh}_4$ , respectively. Figure 4(f) shows the zoomed-in section of the cell marked with a red square in $\text {Mesh}_4$ , which contains an isolated refined cell. The coordinates of the isolated refined cell’s centre $(X_1, Y_1)$ are used to estimate the centre of the vortices. The centre is marked with a blue dot in Figure 4(f) for illustration purposes. The marked red boxes in Figure 4(a)–4(e) demonstrate that the AMR method can capture the location of the centre of vortices within the fourth refined cells.
Erturk and Gokcol [Reference Erturk and Gokcol8] illustrated that the coordinates of centres of location of the primary vortices are stated from the centre of the square cylinder, as shown in Figure 5. Hence, all estimated coordinates for the centres of vortices, $(x_1, y_1)$ , are the results of $(X_1, Y_1)-(10.5, 4)$ , where $(X_1, Y_1)$ is the coordinate of the centre of the isolated refined cell. The estimated centre of vortices for the flow with $Re=5$ is ${(x_1, y_1)=(0.609, 0.172)}$ . We note that the estimated coordinate of the centre of vortices obtained with $\text {Mesh}_4$ is in good agreement with $(x_{\text {ref}}, y_{\text {ref}})=(0.610, 0.180)$ as reported by Erturk and Gokcol [Reference Erturk and Gokcol8].
For the flow with the remaining cases, $10\le Re\le 50$ , a similar mesh refinements pattern, as for the flow with $Re=5$ , was observed, for example, $\text {Mesh}_4$ was refined mostly near the channel’s top and bottom boundaries, in the vicinity of the square cylinder and the horizontal centre line. It is observed that the x and y values of the coordinates of the centre of vortices gradually increase as the Reynolds numbers increase. The setting for the flows in this study is symmetric, but the refined meshes are not. For all the cases, the lack of symmetry in the refined meshes indicates the nonsymmetrical profile of the velocity field. Hence, the calculated velocity fields are not accurate enough. The number of nodes in $\text {Mesh}_4$ for the flow with $10\le Re\le 50$ ranged from 12 989 to 13 914, which were relatively similar to the case with $Re=5$ .
3.2. Estimated location of centre of vortices
Table 1 presents the number of nodes in the fourth refined meshes and the estimated coordinates of vortex centres for the eight test cases. The table also compares the estimated coordinates with the reference centre locations, as reported by Erturk and Gokcol [Reference Erturk and Gokcol8]. For each case, the errors $e_x=x_1-x_{\text {ref}}$ and $e_y=x_1-y_{\text {ref}}$ were found to be within 5%. The errors relative to the reference vortex centres, computed as
were less than 2.5% for all cases. The results obtained by the present method are in good agreement with those obtained by Erturk and Gokcol [Reference Erturk and Gokcol8] using a finer mesh. Moreover, we note that further mesh refinements may be needed to obtain coordinates as close as the actual locations of the vortex centres. Figure 6 shows streamline contours and the flow features near the vortex centre for the case with $Re=5$ and $Re=50$ . The red dots in Figure 6 indicate the location of the centre of vortices.
4. Conclusion
In this study, we verify the accuracy of the 2D AMR with four refinements using the steady flow past a square cylinder. The accuracy is demonstrated by the estimated location of the centre of vortices contained in the first, second, third and fourth refined cells. Since the 2D AMR can be applied a finite number of times, we can perform more refinements if the results are not accurate enough for some cases.