Abstract A finite-element model based on the penalty function formulation of the Navier-Stokes equations governing flows of unsteady, incompressible, non-isothermal, non-Newtonian fluids in three-dimensional enclosures is presented. Power-law and Carreau constitutive relations are used, and the viscosity is assumed to be temperature dependent. The resulting non-linear equations are solved by Picard's method (i.e. direct iteration). The finite-element model is used to analyze several problems, such as flows through circular pipe, cubical cavity, and tubular and square contractions. Wherever possible, comparisons are made between present numerical solutions and available analytical or numerical solutions. The present results are in good agreement with other solutions. The finite-element model developed here can be modified to accommodate other forms of constitutive relations.