**1. Introduction**

The purpose of airfoils is to produce lift with minimum drag. To obtain this result, aerodynamic sections are carefully designed. Very common, in a wide range of fields such as airplanes, ships, sports cars, and sailboats, are the NACA sections [1]. Over the years, many studies have been conducted of these profiles to optimize their shape and also to guide the designers in their choice of the best profile for a certain application. Parallel-sided sections are a particular family of airfoils where the main part of the section has parallel sides, while only the front and aft parts are streamlined. They find applications for sailing dinghies, such as 470, 420 and Optimist, where class rules [2–4] define how the centerboard must be built. Other applications are supporting struts in a fluid flow, and a wind-tunnel application was presented by Pollock [5]. This is the only investigation on parallel-sided sections known to the authors, and the objective of the present work is to fill this knowledge gap, with an emphasis on sailing boats. For this reason, the parameter ranges are taken from the class rules of the dinghies mentioned. The investigation includes two different conditions: upwind (lift coefficient equal to 0.4) and downwind (angle of attack equal to 0.0◦). The upwind case is of interest for all three classes, while the downwind case is of less interest for the 470 and the 420, since they hoist the centerboard while sailing downwind in most wind speeds. For the Optimist, the centreboard is partly hoisted. However, since this type of section can be found in other fields than sailing, the authors are interested in showing the results of both cases.

By far, the most common tools for designing airfoils over the years have been potential flow/boundary layer methods. Even today, the most widely used software for this purpose: XFOIL [6], is based on these theories. However, airfoil design depends crucially on the prediction of separation, both at the leading and trailing edges, and this is a weak point of boundary layer theory. Furthermore, the viscous/inviscid interaction is difficult to model accurately for thick and separated boundary layers. To increase accuracy, less approximate methods based on the Navier-Stokes equations are preferable, and the only realistic possibility, if thousands of computations are to be carried out, is the Reynolds-averaged Navier–Stokes (RANS) method.

Another critical feature of the numerical method is the transition prediction. By designing sections with its pressure minimum far aft, the laminar boundary layer may be maintained over the main part of the foil, reducing the frictional resistance considerably. This is the feature of the so-called laminar sections in the 6-series of the NACA sections [1], but transition is not only important for the friction. There is a crucial dependence of separation on transition, as laminar boundary layers separate more easily than turbulent ones, where the turbulent mixing transports external flow momentum to the inner parts of the boundary layer. This is important, both for separation bubbles at the leading edge, and for the open separation further back on the foil. The latter effect will be seen very clearly in the results below. For a thorough discussion on the effects of transition and separation on airfoil performance, see Larsson et al. [7].

In order to predict these critical features, attention should be paid to both the numerical and physical modelling of the flow [8]. The former should be addressed through a numerical verification study. In the present work, this is performed through a systematic grid variation. Since the method used is based on the RANS equations, a turbulence model is needed, and the choice may have a significant influence on the separation prediction. The motivation for the present choice is given below. Reasons for the choice of transition model are also given.

The numerical method with its parameter settings is described in Section 2. In Section 3, the verification study is presented, including a systematic grid variation and a numerical uncertainty assessment. Section 4 deals with the validation of the method against a NACA 64-006 section, for which measured data are available in [1], and in Section 5 the systematic investigation is presented. Results are shown in Section 6, and the final conclusions are drawn in Section 7.

#### **2. Numerical Method**

#### *2.1. Equations*

In the present project, STAR CCM+ [9] is used for the CFD computations. The Reynolds-averaged Navier–Stokes equations are solved together with the k-ω SST turbulence model [10]. This is by far the most popular turbulence model for hydrodynamic computations, see e.g., Larsson et al. [11] and Hino et al. [12] It is based on the original k-ω model developed for boundary layers in adverse pressure gradients by Wilcox [13], but improved by Menter [10] to avoid the unphysical influence of the external turbulence. The k-ω SST model is recommended for computations of the present kind in the STAR CCM+ Manual [9].

In STAR CCM+, two transition models are available: the two equation Gamma ReTheta model and the one equation Gamma model. The former has been widely used in aerodynamic and turbomachinery applications [14], but has a few drawbacks, and therefore the more modern Gamma model was developed [15]. Since the drawbacks are eliminated in this model, it was selected for the present investigation.

### *2.2. Boundaries*

Figure 1 shows the computational domain and the related boundary conditions applied for the validation case (NACA64-006) and for the systematic investigation of the parallel-sided sections (the zero angle of attack and the constant lift).

**Figure 1.** Boundary definition.

For the NACA 64-006 section, the ambient turbulence intensity, defined as a fraction of the free stream velocity, is 0.006, taken from the wind tunnel specification in [16]. For the systematic investigation, an estimation of the turbulence level in the sea, presented in [17], is adopted. This value is 0.0025. The turbulence viscosity ratio is the ratio between the turbulent and laminar viscosities, and determines the turbulence length scale. Note that, to maintain the ambient turbulence level throughout the domain, source terms are added to the kinetic energy equation; see [9]. This is an important possibility in STAR CCM+. Tests were done without these source terms, and numerical diffusion caused the turbulence to die out rapidly behind the inlet. There was then zero turbulence level at the position of the foil, and this caused the solution to become more unstable. The inner boundary condition is no-slip.

#### *2.3. Discretization*

For the discretization, the finite volume discretization method is used with a segregated flow solver (SIMPLE algorithm). All computations presented in this paper are carried out in steady mode.

#### *2.4. Grid Settings*

The same grid settings are used for the test case and the systematic investigations. Polygonal and prism layer meshes, with the number of prism layers equal to 100 and prism layer stretching equal to 1.0678, are selected. y+ is set to 0.125. The mesh is generated directly in 2D. In Figure 2, the grid close to the test foil, with close-ups at the leading and trailing edges, is presented. With the present transition model, it turned out that a fluctuating laminar separation bubble occurred on the suction side of the nose (lift case) for all tested grids with 50,000 to 400,000 cells. This bubble disappeared for all finer grids, and is therefore assumed to be unphysical. The transition model thus poses strong requirements on the grid density; 400,000 cells is very much for a 2D case. The presented grid settings were the best configuration obtained after the test case analyses, and for this reason, they were also adopted for the systematic investigation.

**Figure 2.** View of the mesh surrounding the test case section, and close-ups of leading and trailing edges. For the test case, the angle of attack is achieved by rotation of the foil, while for the systematic computations the direction of the inflow is varied to obtain the requested lift coefficient (0.4).

#### **3. Verification**

Once the coarsest grid for a stable solution is identified, the verification process to identify the numerical uncertainty is carried out. Two methods are recommended by the ITTC [18]: the factor of safety method by Xing and Stern [19], and the least square root (LSR) method by Eça and Hoekstra [20]. The advantage of the LSR method is that it takes numerical scatter into account by considering solutions from more than three grids. For more details, see [20]. Therefore, this method is selected. Note that only grid uncertainty is considered. Since iterative convergence is achieved by reducing the residuals by 4–5 orders of magnitude, the iterative error can be neglected

Four different grids, the largest one with 6.5 million cells, are created to perform the verification process with the least square root method [20]. Table 1 and Figure 3 show the results. Hi is the linear cell size ratio between a specific grid (i) and the finest one, U is the numerical uncertainty, and U (%CL) is the uncertainty in percent of the coefficient.


**Table 1.** Numerical uncertainty of lift and drag coefficient for the test case at different grid densities.

The very high number of cells for this 2D study leads to very low levels of uncertainty for the finer grids. However, the coarsest grid, with 0.4 million cells, is chosen for further computations. This decision is taken considering the number of simulations that must be performed during the systematic investigation. Although the numerical uncertainty is larger for this grid (1.4% and 3.5% respectively for lift and drag), the difference between the grids in terms of drag coefficients is very small, as can be seen in Figure 4, where the coarse mesh (0.4 million cells) is compared with a finer

one (1.6 million cells). In this comparison, a family of parallel-sided sections with small thickness (T1), long leading edge (LE2) and long trailing edge (TE2) is taken as an example (see Section 5 for the nomenclature). The angle of attack is equal to zero, and five configurations with the same nose radius, but varying trailing edge angles, are evaluated for both the coarse and the fine grid.

**Figure 3.** Regression curve of lift coefficient (**a**) and drag coefficient (**b**).

**Figure 4.** Comparison of drag coefficients for the coarse grid and a finer one. The independent variable is the trailing edge angle (see below).

#### **4. Validation Case**

To validate the numerical model, an effort was made to find a suitable test case in the literature. The aim was to find a thin section at a relevant Reynolds number. However, although there are several such cases reported, even with detailed flow data, none appeared to present accurate force measurements. Therefore, a NACA 64-006 section is selected as a test case [1]. The lowest Reynolds number in the NACA tests is 3,000,000, while in the present systematic computations the Reynolds number varies from 300,000 to 1,500,000. The results from the validation simulation are compared with the experimental data for angles of attack between 0.0◦ and 4.0◦ in Figure 5 and Table 2. Apparently, the transition model is not accurate enough to predict the width of the drag bucket, but since this study is carried out for 0◦ and approximately 4◦, the less perfect match of results between 1.5◦ and 2.0◦ should not influence the systematic investigation.

**Figure 5.** Comparison between numerical and experimental results. Reynolds number 3,000,000.


**Table 2.** Numerical and experimental results for different angle of attack.

It should be stressed that this is not a formal validation (see, for instance, [18]). For this to be possible the experimental uncertainty is required, and this is not available in [1]. Therefore, this is merely a comparison with experimental data.

#### **5. Systematic Investigation**

The authors are mainly interested in 470, 420 and Optimist centerboards, so the variable ranges are chosen according to the rules of each class [2–4] and their respective typical sailing conditions. However, these values are reasonable for other applications of these sections. Table 3 shows the limits for the Reynolds number (based on chord length and typical speed for each boat), thickness, leading edge length and trailing edge length.


**Table 3.** Ranges of Reynolds numbers, thicknesses, leading edge lengths, trailing edges lengths of 470, 420 and Optimist dinghies.

Table 4 summarizes the setup of the systematic investigation, to cover as closely as possible all cases for the three sailboats, and considering the time to perform the entire investigation. Three Reynolds numbers are chosen, as well as two sailing conditions, two thicknesses, two leading edge lengths, two trailing edge lengths, five nose radii and five trailing edge angles. The latter are defined in terms of a reference angle *AREF*, where

$$A\_{REF} = \arctan\left(1.5\frac{\text{T}}{\text{TE}}\right) \tag{1}$$

**Table 4.** Values and codes for the systematic investigation of Parallel-Sided sections.


For more details, see Appendix A. With all these variables, a set of 200 shapes is generated, and in total, 1200 cases are investigated, considering the combination of the three Reynolds numbers and the two sailing conditions.

For the upwind sailing condition, the angle of attack is not defined. A target value of lift (CL = 0.4) is set instead, so the angle of attack is slightly different between the cases. It is the given side force from the sails that should be balanced, and this will be achieved at different leeway angles for the different sections. The authors choose a lift coefficient equal to 0.4 as a target value after conducting a study of typical leeway angles of the 470 and Optimist dinghies using a velocity prediction program (VPP). For 420 and 470, the 0◦ case may not be relevant, since the centerboard is hoisted while sailing downwind. For the Optimist, the centerboard is usually only partly hoisted during the downwind leg, so the 0◦ case can be useful for this type of dinghy. Furthermore, the 0◦ case for this type of sections is of interest in other fields than sailing, where the angle of attack is often 0◦.

#### **6. Results**

In the following figures, the computed results are presented. To enhance readability, only the best case for each family of thickness (T), leading edge length (LE) and trailing edge length (TE) is represented, since the trend is the same inside each family, and the results are very close. In particular, for 0◦ angle of attack, the nose radius has a very small effect on the drag, and the same is true for the trailing edge angle at constant lift. By presenting the results in this way, the authors believe that it is

easier to compare the results and draw conclusions. The drag coefficient for the cases with a Reynolds number of 300,000 are averaged over the last 1000 iterations, due to a slightly oscillating solution.

First, the results of all families are presented. Figures 6–8 show the results for the Reynolds numbers 300,000, 900,000 and 1,500,000, respectively. The top left plot in each figure is for the thin sections (4%) at 0◦ angle of attack, while the top right figure shows the corresponding results for the thick sections (8%). The two lower figures show the results for the CL = 0.4 case.

**Figure 6.** Drag coefficient (CD <sup>×</sup> 103) for Reynolds number 300,000. Independent variable: trailing edge angle (**a**,**b**); nose radius (**c**,**d**).

**Figure 7.** Drag coefficient (CD <sup>×</sup> 103) for Reynolds number 900,000. Independent variable: trailing edge angle (**a**,**b**); nose radius (**c**,**d**).

**Figure 8.** Drag coefficient (CD <sup>×</sup> 103) for Reynolds number 1,500,000. Independent variable: trailing edge angle (**a**,**b**); nose radius (**c**,**d**).

Having presented the results for all families of sections, the results will now be analyzed and commented on. In Section 6.1 results for zero angle of attack are discussed, while in Section 6.2 the results for CL = 0.4 are considered. In order to explore the physics, field plots of the axial velocity are presented, together with line plots of the pressure or friction variation along the section on the two sides. To identify regions of separation, the lower limit of the axial velocity is set to zero, such that negative velocities, i.e., separations, are indicated by a white color.

#### *6.1. Zero Angle of Attack*

It is obvious from Figures 6–8 that the thinner sections have considerably lower drag than the thick ones at zero angle of attack. A comparison is made in Figure 9, where the thicker profile exhibits separation (white area) at the trailing edge. This causes a reduced pressure in the build-up region between X/L from 0.6 to 0.9. The reduced pressure in this region increases drag, although the pressure is higher at the trailing edge. T2 has a lower pressure also further forward, but this is mostly in the parallel part and has a small effect on the pressure drag.

**Figure 9.** Plots of axial velocity (upper) and pressure distribution (lower) for a thin and a thick section.

As it appears from Figures 6–8, the length of the leading edge is unimportant for thin sections. For thick sections, it is important in some cases. Figure 10 shows that the profile with a short leading edge performs better than the one with the longer leading edge. This is due to the different nature of the boundary layers. As can be seen in the figure, the skin friction coefficient for the short leading edge is much larger than for the other case. This indicates that the flow on the short leading edge section is turbulent, while it is laminar on the other section. The reason for this is that the pressure minimum at the forward shoulder is much lower for the short leading edge, causing a large positive pressure gradient behind the minimum, thus promoting both transition and separation. Although the laminar flow gives less friction, the drag of the section is high, since a larger (unsteady) separation zone is created at the trailing edge of this foil. A laminar boundary layer separates much more easily than a turbulent one, where energy is transferred through the boundary layer due to the stirring effect of the turbulence. Note that the separated region on the short leading edge foil between X/L = 0.1 and 0.2 has negative friction, indicating a small separation bubble in this region. Since this separation is on the parallel part, it does not increase drag. The long leading edge has an unsteady separation, so the figures show instantaneous results.

Figures 6–8 show that the trailing edge should be as long as possible. As can be seen in Figure 11, the short trailing edge has massive separation, while no separation can be seen on the long trailing edge. The separation reduces the pressure from about X/L = 0.7 and backwards, which results in larger drag.

**Figure 10.** Velocity plot and skin friction coefficient of short leading edge length (**left**) and long leading edge length (**right**).

**Figure 11.** Velocity plot and pressure plot of short trailing edge (**left**) and long trailing edge (**right**) sections.

The nose radius has a very small effect on the drag. An example is given in Figure 12, where the drag for all thin sections with long leading and trailing edges are shown for the intermediate Reynolds number. Plotting drag versus trailing edge angle, it can be seen that the curves for all five nose radii collapse.

The trailing edge angle should be small for the thin sections. For the thick sections, this is also true, however with one exception: the family with a short trailing edge. In this case, the optimum solution is not obtained with the smallest trailing edge angle, but with a medium angle. Figure 13 shows results for three trailing edge angles. The flow for all three cases is laminar over the entire profile, so when it gets close to the trailing edge, it separates. When the angle is large (figure to the left), there is a large curvature at the aft shoulder which promotes (unsteady) separation, and when the angle is small (figure on the right), the aft end exhibits an inflexion, which causes a large curvature slightly more forward. This again promotes separation. The smallest curvature is obtained for an intermediate angle (middle figure), which has a slightly smaller separated area and a smaller drag. Note that the separation is unsteady, so the plots show results at one time instant.

**Figure 12.** Five different radii are compared. *x*-axis trailing edge angle, *y*-axis drag coefficient times 103.

**Figure 13.** Velocity plot and skin friction coefficient of large (**left**), medium (**center**) and small (**right**) trailing edge angle sections.

The best thin parallel-sided sections have on average 13% higher drag than the four-digit NACA section with the same thickness ratio. For the thick sections, the average drag increase is 30%.
