# Numerical Simulation of KBKZ Integral Constitutive Equations in Hierarchical Grids
## 4 Results
### Numerical Simulation of 4:1 Abrupt Planar Contraction Problem  
In this section, we show the simulations for the 4:1 abrupt planar contraction flow. This problem is interesting because, for instance, the flow near the contraction is a complex mixture of shear and elongation, and secondary fluid motions might exist, even in the Newtonian limit (see [15]). For this reason, contraction flows have been extensively studied previously in the literature (see [3; 14; 15; 26]). Figure 8 shows the domain representation, where we adopted a dimensionless parabolic inlet velocity profile \(u(y)=\frac{3}{8}\frac{1}{4}(2-y)(2+y)\). The scaling parameter \(L=0.0064\) m is the height of the small channel, and we used \(N=50\) deformation fields, \(\epsilon=0.1\) and the time interval \(s_{c}=0.1\) s. In Table 6, we report the scaling parameters of average velocity \(\bar{u}\) used in all simulations and the dimensionless parameters \(Re=\frac{\rho LI}{\eta_{0}}\), \(De=\frac{\lambda LI}{L}\), \(De(\dot{\gamma})\) (see [15]) and the characteristic shear rate \(\dot{\gamma}\). The mesh \(M_{1}\) used in these simulations is shown in Figure 9. We used three levels of refinement, with the most refined part near the contraction region. We also used one uniform mesh \(M_{2}\) for \(De=0.94\) in order to check the convergence solutions in two meshes. In \(M_{2}\), \(d_{x}=d_{y}=0.03125\) m, and in \(M_{1}\), we use small values of \(dx=dy=0.03125\) m as well as larger values, \(dx=dy=0.125\) m. \begin{table}
\begin{tabular}{c c c c c} \hline \hline \multicolumn{5}{c}{**Planar Contraction Flows**} \\ \hline \(\bar{u}[\frac{m}{s}]\) & \(\dot{\gamma}[\frac{1}{s}]\) & \(De(\dot{\gamma})\) & \(Re\) & \(De\) \\ \hline
0.044 & 13.9 & 0.38 & 0.16 & 0.41 \\
0.100 & 31.3 & 0.55 & 0.36 & 0.94 \\
0.150 & 48.4 & 0.66 & 0.56 & 1.45 \\
0.221 & 69.1 & 0.77 & 0.80 & 2.07 \\ \hline \hline \end{tabular}
\end{table} 
Table 6: Flow parameter values used in the contraction problem. Figure 8: Domain representation. Figure 9: Graphical representation of the mesh \(M_{1}\) used in the contraction simulation. In Figure 10, we illustrate the centerline axial profile velocity \(u_{x}(y)\) solutions for different values of the Deborah number \(De=0.41,0.94,1.45\) and \(2.07\) for the HiGFlow (green lines) and OpenFOAM (black lines) systems using mesh \(M1\). The experimental results of [PERSON] et al. [15] and numerical results of [PERSON] [14] and [PERSON] et al. [3] (for \(De=1.45\)) are also shown for comparison purposes. For the case with \(De=0.94\), we show two solutions using our methodology in mesh \(M_{1}\) (non-uniform mesh) and mesh \(M_{2}\), where we can see good agreement between the solutions. For this reason, we adopted the \(M_{1}\) mesh to simulate the other cases with different values of number \(De\). For the profile \(u_{x}(y)\), the HiGFlow system showed good agreement with the OpenFOAM solutions in the regions before and after the contraction. Near to contraction (see Figure 10), we have a region of instability and the methodology behaves differently, but our results showed similar behavior to the instabilities presented in the works of [PERSON] et al. (orange triangles) [15], [PERSON] [14] (blue bullet) and [PERSON] et al. (violet square) [3]. In Figure 11, we show the numerical solution using HiGFlow (green lines) and OpenFOAM (black lines) systems for the tensor components \(\tau_{xx}\) and \(\tau_{yy}\) with the same flow parameter values reported in Table 6 using \(M_{1}\) and \(M_{2}\) for the case with \(De=0.94\). The methodologies presented in this work are different. OpenFOAM uses the finite volume method while HiGFlow approximates the equations using finite differences. Therefore, the solutions obtained will not be equal but should be comparable. Outside the contraction region, the OpenFOAM and HiGFlow solutions are very similar for all values of \(De\).