Magnetoelectric materials have received astonishing attention due to their great potential for nanotechnology and magnetoelectric applications. Herein, the ferroelectric, magnetic and magnetoelectric properties of multiferroic BiFeO3 (BFO) thin films have been carried out using Monte Carlo simulation investigations in the framework of the Heisenberg model. Exchange coupling interactions in magnetic and ferroelectric sublattices that correspond to the experimental critical temperature were estimated. Temperature dependence of the internal energy, specific heat, magnetization, electric polarization and their susceptibilities in BFO thin films have been systematically analyzed. Moreover, the effect of magnetoelectric coupling interaction Jme on M-H and P-H hysteresis loops and magnetoelectric voltage coefficient were studied. A large magnetoelectric voltage coefficient of 104 mV/(cm.Oe) was predicted in BFO thin films with 25% defects. Furthermore, the defect mechanism can be used to control the polarization switching by the magnetic field in BFO system. The obtained results highlight the large magnetoelectric effect in BFO multiferroic thin films for multifunctional magnetoelectric devices. \textcopyright 2021 Elsevier B.V.Read less <