Nonlinear relationship between Reynolds stresses and the rate of strain of nonlinear k-$varepsilon$models is evaluated theoretically by using the boundary layer assumptions against the turbulence-driven secondary flows in noncircular ducts and then their prediction performance is validated numerically through the application to the fully developed turbulent flow in a square duct. Typical predicted quantities such as mean axial and secondary velocities, turbulent kinetic energy and Reynolds stresses are compared with available experimental data. The nonlinear k-$varepsilon$ model adopted in a commercial code is found to be unable to predict accurately duct flows with the prediction level of secondary flows one order less than that of the experiment.