FSI simulations of flapping motions have been widely investigated to develop a flapping-wing micro air vehicle. Because an intensive parametric study is important for the product design, a computationally efficient model is required. The purpose of the present study was to develop a reduced-order model of flapping motion. Among the various methods available to solve FSI problems, we employed the Dirichlet–Neumann partitioned iterative method, in which three sub-systems (fluid mesh update, fluid analysis, and structural analysis) are executed. In the proposed analysis system, first, snapshot data of structural displacement, fluid velocity, fluid pressure, and displacement for the fluid mesh update were collected from a high-fidelity FSI analysis. Then, the snapshot data were used to create low-dimensional surrogate systems of the above three sub-systems based on the POD under Galerkin projection (i.e., the POD-Galerkin method). In numerical examples, we considered a two-dimensional FSI problem of simplified flapping motion. The problem was described via two parameters: frequency and amplitude of flapping motion. We demonstrated the effectiveness of the presented reduced-order model in significantly reducing computational time while preserving the desired accuracy.