Also the inverse transformation can be done using the following relation
From these relations, we can calculate the transformation metrics as the following
The 2-D steady incompressible flow inside a triangular cavity is governed by the Navier-Stokes equations. We consider the N-S equations in streamfunction and vorticity formulation, such that
where y is streamfunction and w is vorticity, x and y are the Cartesian coordinates and Re is the Reynolds number. We note that these equations are non dimensional and a length scale of (h - cy)/3 and a velocity scale of U, ie. the velocity of the lid, are used to non-dimensionalize the parameters and the Reynolds number is defined accordingly.
Using the chain rule, the governing N-S equations in general curvilinear coordinates are the following
Substituting for the transformation metrics in Equations (3 and 4), we obtain the equations that govern the flow in a triangular cavity shown in Figure 1 as the following
The boundary conditions in x- and y-plane are the following
where n is the unit normal vector and t is the tangent vector in clockwise direction and (u,v) are the velocity components in x- and y-direction where
Therefore, in the computational domain, the boundary conditions on the side A-B are
On sides A-C and B-C we have
Also note that, on side B-C one can also show that the derivative of the streamfunction in the wall normal direction would be homogenous.
Erturk et al.  stated that for square driven cavity when fine grids are used, it is possible to obtain numerical solutions at high Reynolds numbers. At high Reynolds numbers thin boundary layers are developed along the solid walls and it becomes essential to use fine grid meshes. Also when fine grids are used, the cell Reynolds number or so called the Peclet number, defined as Rec = uDh/n, decreases and this improves the numerical stability (see Weinan and Jian-Guo  and Tennehill et al. ). In this study, we use a fine grid mesh with (512×512)/2 grid points. On this mesh, we solve the governing equations (10 and 11) using SOR method. For the vorticity values at the boundary points we use Thom's formula . We note an important fact that, it is well understood (see Weinan and Jian-Guo , Spotz , Huang and Wetton , Napolitano et al. ) that, even though Thom's method is locally first order accurate, the global solution obtained using Thom's method preserve second order accuracy O (Dh2). On side A-B the vorticity is equal to
where the subscripts i and j are the grid indices in x- and h-directions respectively and the grid index 0 refers to points on the wall and 1 refers to points adjacent to the wall and also Dh refers to grid spacing.
Similarly on side A-C the vorticity is equal to
and on side B-C, it is equal to
3 RESULTS AND DISCUSSIONS
During the iterations, as a measure of convergence, we monitored several residuals. We define the residual of the steady state streamfunction and vorticity Equations (10 and 11), as the following
These residuals indicate the degree to which the numerical solution has converged to steady state. In our computations, in all of the triangle geometries considered, at any Reynolds numbers, we considered that the convergence was achieved when the absolute values of the maximum of the residuals defined above in the computational domain ( max( |Ry| ) and max( |Rw| ) ) were both below 10-10. Such a low convergence level is chosen to ensure the accuracy of the solution. At these convergence levels, the maximum absolute change in the streamfunction variable in the computational domain ( max( |yi,jn + 1 - yi,jn| ) ) was in the order of 10-16 and for the vorticity variable ( max( |wi,jn + 1 - wi,jn| ) ) it was in the order of 10-14. In other words, our numerical solutions of the streamfunction and the vorticity variables are accurate up to 15 and 13 digits respectively. Also at these convergence levels, the maximum absolute normalized change in the streamfunction and vorticity variable in the computational domain ( max( |(yi,jn + 1 - yi,jn )/ yi,jn| ) and max( |(wi,jn + 1 - wi,jn ) / wi,jn| ) ) were in the order of 10-13 and 10-12 respectively. In other words, at convergence the streamfunction variable was changing with 10-11 percent of its value and the vorticity variable was changing 10-10 percent of its value, at each iteration step. These numbers ensure that our numerical solutions are indeed very accurate.
In this study we considered different triangle geometries. In order to demonstrate how we define the Reynolds number, let us consider a dimensional equilateral triangle with coordinates of corner points
If we use the top wall velocity U and the length scale a in nondimensionalization, then the Reynolds number is defined as Ua/n and the coordinates of the corner points of the nondimensional equilateral triangle become
We, first considered this equilateral triangle geometry which was also considered by McQuain et al. , Ribbens et al. , Jyotsna and Vanka , Li and Tang  and Gaskell et al. . We solved the flow in this triangle cavity at various Reynolds numbers ranging between 1 and 1750. We note that, if the length of one side of the triangle () was used in nondimensionalization, as it was used by Gaskell et al. , then our Reynolds number of 1750 would be fold such that it would correspond to a Reynolds number of 6062. Qualitatively, Figure 2 shows the streamline contours of the flow at various Reynolds numbers and also Figure 3 shows the vorticity contours at the same Reynolds numbers.
In terms of quantitative analysis, Table 1 tabulates the center of the primary eddy and the streamfunction and vorticity values at the core, together with results found in the literature. Our results agree well with McQuain et al. , Ribbens et al.  and agree excellent with Gaskell et al.  up to the maximum Reynolds number (Re=500) they have considered. However, the results of Li and Tang  differ from our results.
The difference can be seen more obvious in Figure 4, where we plot the vorticity values at the center of the primary vortex tabulated in Table 1, with respect to the Reynolds number. The results of Li and Tang  start to behave differently starting from Re=100 from the rest of the results.
Also in Figure 4, the results of McQuain et al.  and Ribbens et al.  start to deviate from present results and Gaskell et al.'s  results after Re=200. Having a different behavior, the vorticity value of McQuain et al.  and Ribbens et al.  shows an increase such that their vorticity value at Re=500 is greater than the vorticity value at Re=350, whereas the present results and Gaskell et al.'s  results show a continuous decrease. We believe that these behaviors are due to the coarse grids used in those studies and in order to resolve these behaviors we decided to solve the same flow problem using several coarse grid meshes. We computed the flow using (64×64)/2, (128×128)/2 and (256×256)/2 grid points also. For a given grid mesh we have solved the flow for increasing Reynolds number until a particular Re where the solution was not convergent but oscillating. For this particular Re when the number of grid points was increased, the convergence was recovered and we were able to obtain a solution.
In order to demonstrate the effect of the grid mesh on the numerical solution, in Figure 5 we plotted the vorticity value at the center of the primary vortex obtained using different grid mesh, with respect to Reynolds number.
From the results in Figure 5, we conclude that (64×64)/2 grid mesh is too coarse for this flow problem. In this figure, looking at the solution of (128×128)/2 and (256×256)/2, the vorticity value decreases with increasing Reynolds number however shows an increase at the largest Reynolds number where beyond that Re the solution is not converging. We note that this behavior is similar to that of McQuain et al.  and Ribbens et al.  in Figure 4.
In Figure 5 looking at the solutions of different grid mesh, for example at Re=350, we see that as the number of grid points in the mesh is increased, the vorticity value at the center of the primary eddy is decreased however this decrease slows down as Dhgets smaller. This is not surprising since the solution gets more accurate as the number of grid points is increased, and also since our solution is spatially second order accurate, O(Dh2), the computed vorticity value should be expected to converge to the real value with Dh2. At the finest grid mesh, (512×512)/2, at high Reynolds numbers our the vorticity value seems to stagnate around 1.16. We note that, we expect the vorticity value to be a little lower than this value, if the number of grid points is increased furthermore. However increasing the number of grid number furthermore increases the computing time, therefore we will not consider a finer grid mesh than (512×512)/2. This grid study shows that fine grids are necessary for accurate solutions for triangular cavity flow especially at high Reynolds numbers.
McQuain et al.  analytically obtained the constant vorticity value of the inviscid rotational flow at infinite Reynolds number as =1.054 using Batchelor’s mean-square law. This infinite Reynolds number value of the vorticity of 1.054 is shown with the dotted line in Figure 4. Our computed primary eddy vorticity asymptotes to a value around 1.16 as the Reynolds number increases, which is a little higher than the theoretical value. We note that the theory, assumes that at infinite Reynolds number the whole inviscid fluid rotates as a solid body with a constant vorticity which is coupled to boundary layer flows at the solid surface. However in reality, from Figure 2 it is evident that the whole fluid does not rotate as a solid body, instead there appears progressively smaller counter-rotating recirculating regions at the bottom corner and also one towards the upper left corner. We note that, especially the eddies at the bottom corner occupy some portion of the corner as the Reynolds number increases. However, the increase in the size of the portion of the bottom corner eddies almost stops after Re=500 and the size of the primary eddy remains almost constant beyond and so is the value of the vorticity at the core of the primary eddy. It looks like, for an equilateral triangular cavity flow at high Reynolds numbers, the mean square law predicts the strength of the primary eddy within an error due to the effect of the secondary eddies. For circular or elliptic boundaries Ribbens et al. , for square cavity flow Erturk et al. , and for rectangle cavities McQuain et al.  have shown that the mean square law is approximately valid and successful in predicting the strength of the primary eddy at high Reynolds numbers. Due to small corner angles in a triangle geometry, for triangular cavities the mean square law is not as successful as it was in other geometries.
We then considered an isosceles triangle which was also considered by Jyotsna and Vanka , where
We note that, our definition of Reynolds number is equivalent to one fourth of the Reynolds number definition used by Jyotsna and Vanka .
Figure 6 shows the streamline contours of the flow in this triangle at various Reynolds numbers. Table 2 compares the location of the primary eddy center with that of Jyotsna and Vanka  and Gaskell et al. . The results agree with each other, although we believe that our results are more accurate.
Moffatt  have analytically studied the recirculating eddies near a sharp corner in the Stokes regime and predicted that the distance from the corner to the centers of eddies and also the velocity at the dividing streamline between the eddies, follow a geometric sequence.
Using the considered triangle geometry, Jyotsna and Vanka  have used the small Reynolds number solutions to compare with Moffatt's  Stokes regime predictions. However, when Reynolds number increases the asymmetry in the flow field increases and the comparison with the Stokes regime could not be possible. In order to be able to compare our results with Moffatt , we solve the considered triangle geometry for the Stokes regime, i.e. Re=0. The solution at Re=0, as seen in Figure 6a, is symmetrical. For this Stokes regime solution, we first calculate the u-velocity along a perpendicular line from the bottom corner to the top wall. In this velocity profile, the points where the u-velocity is zero correspond to the eddy centers. The distance of these points to the bottom corner would give us rn, where we use in calculating the ratio rn /rn + 1. In the velocity profile, the points where the u-velocity has local maximum correspond to the u-velocity at the dividing streamline between the eddies, where Moffatt  have used these velocities as a measure of the intensity of consecutive eddies. These local maximum u-velocities would give us In, where we use in calculating the ratio In /In + 1. For the considered triangle geometry where q=28.072°, Table 3 tabulates the calculated ratios of rn /rn + 1 and In /In + 1 along with analytical predictions of Moffatt  and the agreement is good.
We then considered an isosceles triangle which was also considered by Gaskell et al. , with
where q=20°, 30°, 40°, 50°, 60°, 70°, 80° and 90°. We will consider Stokes flow in these triangles, such that we solve for Re=0. The analytical predictions of Moffatt  for these θ values, are tabulated in Table 4.
Figure 7 qualitatively shows the change in the Stokes flow in an isosceles triangle as the corner angle, θ, changes. For this case, our computed relative positions and intensities of eddies are tabulated in Table 5 along with results of Gaskell et al. . The agreement between our computed results and that of Gaskell et al.  and also predictions of Moffatt  is good.
We note that, Moffatt's  analytical analysis predicts that there will be infinite sequence of eddies at the bottom corner. This result is, of course, impossible to capture for a numerical analysis due to the resolution of the grid mesh. However since we have used a finer grid mesh, with (512×512)/2 grid points, than the grid mesh used by Gaskell et al. , with 6605 nodes, we were able to capture more sequence of eddies than Gaskell et al. , at the sharp bottom stagnant corner.
We then considered an isosceles right triangle with the 90° corner being at the top right corner, such that with corner points
The flow in this considered triangle is solved up to Re=7500. Figure 8 shows the flow topology as the Reynolds number increases and Table 6 tabulates the properties of the primary eddy, the streamfunction and the vorticity value at the center of the eddy and also the location of the center, for future references.
We then also considered an isosceles right triangle with the 90° corner being at the top left corner, such that with corner points
Using this triangle geometry, we were able to obtain solutions up to Re=2500. Figure 9 shows the flow topology as a function of Reynolds numbers. In Table 7, for the considered triangle geometry, we tabulated the properties of the primary eddy, i.e. the eddy that is closest to the moving lid, for future references.
As it is obvious from Figures 8 and 9, in both cases the flow behaves very differently as the Reynolds number increases, which shows that the flow structures in a triangle cavity are greatly affected by the triangle geometry.
In this study we have presented high accurate, fine grid solutions of 2-D steady incompressible flow in triangle cavities. The governing equations are solved up to very low residuals at various Reynolds numbers. The computed solutions are compared with numerical solutions of McQuain et al. , Ribbens et al. , Jyotsna and Vanka , Li and Tang  and Gaskell et al.  and also with analytical predictions of Moffatt , and good agreement is found. Our results showed that for an equilateral triangular cavity flow, Batchelor's mean-square law is not as successful as it was in square or rectangle cavity flows, due to small stagnant corner angle.
File translated from TEX by TTH, version 3.63.