AbstractA finite element model is developed based on the penalty formulation to study incompressible laminar flows. The study includes a number of new quadrilateral and triangular elements for 2‐dimensional flows and a number of new hexahedral and tetrahedral elements for 3‐dimensional flows. All elements employ continuous velocity approximations and discontinuous pressure approximations respecting the LBB condition of numerical instability. An incremental Newton–Raphson method coupled with the Broyden method is used to solve the non‐linear equations. Several numerical examples (colliding flow, cavity flow, etc.) are presented to assess the efficiency of e