Recent observations indicate that hydrogen-poor superluminous supernovae (SLSNe-I) often display bumpy declining light curves. However, the cause of these undulations remains unclear. In this paper, we have improved the magnetar model, which includes flare activities. We present a systematic analysis of a well-observed SLSN-I sample with bumpy light curves in their late phase. These SLSNe-I were identified from multiple transient surveys, such as the Pan-STARRS1 Medium Deep Survey and the Zwicky Transient Facility. Our study provides a set of magnetar-powered model light curve fits for five SLSNe-I, which accurately reproduce the observed light curves using reasonable physical parameters. By extracting the essential characteristics of both the explosions and central engines, these fits provide valuable insights into investigating their potential association with gamma-ray burst engines. We found that the SLSN flares tend to be the dim and long extension of GRB flares in the peak luminosity versus peak time plane. Conducting large-scale, high-cadence surveys in the near future could enhance our comprehension of both SLSN undulating light curve properties and their potential relationship with GRBs by modeling their light curve characteristics.