## 2 Numerical Method
\[-\omega=\frac{\partial^{2}\Psi}{\partial X^{2}}+\frac{\partial^{2}\Psi}{ \partial Y^{2}} \tag{1}\]  
\[U\frac{\partial\omega}{\partial X}+V\frac{\partial\omega}{\partial Y}=-\frac{ \varepsilon^{2}}{\mathrm{Re}\cdot Da}\omega-\frac{\varepsilon^{2}C_{F}U_{M}}{ \sqrt  
The relevant empirical coefficients in the present porous medium of pin-fins, such as \(K\), \(C_{F}\), \(k_{s}\)\({}^{*}\), \(k^{*}\) and \(h_{fs}\), generally do not have universal values. This work lists the porous properties of pin-fin heat sinks studied herein in Table 1 by using the relevant empirical formula reported in the open literatures [18, 22, 23]. Some heat sinks in the present study are the same as those used in previous study [18]. This is for verification purposes since the single heat sink in the previous study [18] might be like the first row of heat sinks in the present study. The \(h_{fs}\) of the square pin-fin heat sink with uniform in-line arrangement is predicted by employing the empirical equations suggested by [PERSON] _et al._[24]:  
\[h_{fs}= \frac{\varepsilon k_{f}}{d}\left(0.36283\,\mathrm{Re}_{d}^{\ 0.54219} \right)\ \mathrm{for}\ Re_{d}<1000 \tag{7}\] \[h_{fs}= \frac{\varepsilon k_{f}}{d}\left(0.04433\,\mathrm{Re}_{d}^{\ 0.82934} \right)\ \mathrm{for}\ Re_{d}\geq 1000 \tag{8}\]  
where \(Re_{d}\) is defined as \(\rho u_{max}d/\mu\); \(u_{max}\) is the average maximum velocity between pin-fins, and \(d\) is the pin-fin thickness. Figure 3a schematically depicts the numerical domain and the boundary conditions. The relevant assumptions of boundary conditions for the vorticity (\(\omega\)) and the stream function (\(\Psi\)) refer [PERSON] [25], who brilliantly discussed various boundary conditions for \(\omega\) and \(\Psi\). The fluid is assumed to enter the inlet with a uniform temperature (\(T_{i}\)) and a uniform velocity (\(u_{i}\)). The top and bottom faces are symmetric to the adjacent zones (not included in the present computational domain). The fluid phase at the outlet face is assumed to meet the zero diffusion condition since the downstream length was chosen to be sufficiently long. This work applies the power-law developed by [PERSON] [26] to disperse the equations and employs the SIS (Strong Implicit Solver) algorithm proposed by [PERSON] [27] to solve the related dispersal equations. All of the resolutions undergo the grid independence test and convergence test. The grid systems are separated into various classes to fit the physical model. The grid points, as shown in Figure 3b, are 1224 \(\times\) 49-1440 \(\times\) 77. The grid size falls as \(X_{T}\)\({}^{*}\) and \(X_{L}\)\({}^{*}\) values decrease. The \(X_{T}\)\({}^{*}\)(= \(S_{T}\)\({}^{*}\)/\(L\)) and \(X_{L}\)\({}^{*}\)(= \(S_{L}\)\({}^{*}\)/\(L\)) are the relative center-to-center longitudinal and transverse distances between the adjacent heat sinks, respectively. The iteration ends when the variables meet the criterion,  
\begin{table}
\begin{tabular}{c c c c c c c c} \hline \hline \multicolumn{6}{c}{**In-Line Square Pin-Fin Heat Sink with Uniform Distributions (\(S_{L}\) = \(S_{T}\)), \(L\) = 0.067 m, \(H\) = 0.0237 m**} \\ \hline
**Test Specimens** & **Sample 1** & **Sample 2** & **Sample 3** & **Sample 4** & **Sample 5** & **Sample 6** & **Sample 7** & **Sample 8** \\ \hline \(d\) (m) & 0.00596 & 0.00496 & 0.00372 & 0.00479 & 0.0067 & 0.00298 & 0.00248 & 0.00305 \\ \(n\) & 9 \(\times\) 9 & 9 \(\times\) 9 & 9 \(\times\) 9 & 7 \(\times\) 7 & 5 \(\times\) 5 & 9 \(\times\) 9 & 9 \(\times\) 9 & 11 \(\times\) 11 \\ \(\varepsilon\) & 0.358 &_Appl. Sci._**2015**, \(5\)  
\[\sum_{i,j}\left|\frac{F_{i,j}^{r+1}-F_{i,j}^{r}}{F_{i,j}^{r+1}}\right|\leq 1 \times 10^{-5} \tag{9}\]  
where \(F\) represents \(\Psi\), \(\omega\) or \(\theta\). The subscripts \(i\) and \(j\) refer to the \(i\)-th and the \(j\)-th grid-nodes in the \(X\) and \(Y\) directions, respectively. The superscript \(r\) indicates the \(r\)-th iteration.