<p>This study presents a computationally scalable framework using the lattice Boltzmann method combined with multi direct-forcing immersed boundary scheme and topology-confined block refinement to simulate complex flows past two- and three-dimensional flapping wings at large pitching and flapping frequencies. The framework is validated and analyzed by using Kolmogorov scales and vortex identification methods to correlate instantaneous coherent structures with flow physics. For a two-dimensional flapping plate with high pitching amplitude, the present numerical solver obtains identical leading edge vortex and lift coefficient compared to the experimental results. For the three-dimensional flapping wing, the Liutex identification method is used for the first time, classifying rotational vortex structures into leading edge vortex, tip and root vortex, trailing edge vortex, and U-shaped vortex (combination of the leading edge, tip, and trailing edge vortices). The formation of a U-shaped vortex is significantly correlated with the jump of lift coefficient; where the maximum lift coefficient at the flapping frequency of 0.4 jumps four times higher than that at flapping frequency of 0.2. At low flapping frequency, only the leading edge and tip vortex appear. However, at high flapping frequency, the leading edge vortex detached from upper wing surface because of transverse velocity; while the trailing edge vortex terminated its connection with a U-shaped vortex due to helical motion of the tip vortex. In addition, the leading edge vortex of high vortex strength strongly maintained its connection to the U-shaped vortex; and the small vortex strength of the trailing edge vortex also causes instability of the U-shaped vortex.</p>