A numerical method is developed to solve the two-dimensional incompressible viscous flow with rotating boundaries of arbitrary shape using the boundary-fitted curvilinear coordinate transformation technique. The stream function-vorticity formulation is developed in which the requirement of the single-valuedness of the pressure is explicitly included to determine the boundary values of the stream function. The method is verified by the analysis of the flow between two long, eccentric cylinders. The method is also applied to the transient flow analysis in which the flow is driven by arbitrarily shaped rotating blades, and it is confirmed that the flow patterns can be successfully simulated using the proposed method even if the system has a moving boundary.