In this paper, the coupled fractional Ginzburg-Landau equations are first time investigated numerically. A linearized implicit finite difference scheme is proposed. The scheme involves three time levels, is unconditionally stable and second-order accurate in both time and space variables. The unique solvability, the unconditional stability and optimal pointwise error estimates are obtained by using the energy method and mathematical induction. Moreover, the proposed second-order method can be easily extended into the fourth-order method by using an average finite difference operator for spatial fractional derivatives and Richardson extrapolation for time variable. Finally, numerical results are presented to confirm the theoretical results.The Ginzburg-Landau equation (GLE) has been used to model a wide variety of physical systems [1]. The existence, regularity, and the long-time behavior of the exact solution for the GLE are given in [2][3][4][5]. For the global well-posedness and the existence of the global attractor of the GLE, readers can refer to Refs. [6,7]. Sakaguchi and Malomed [8] proposed the coupled Ginzburg-Landau equations (CGLE) to describe the Bose-Einstein condensates and nonlinear optical waveguides and cavities. And the linearly coupled complex cubic-quintic Ginzburg-Landau equations [9] and the coupled Ginzburg-Landau equations with higher-order nonlinearities [10] are also discussed.The fractional generalization of the GLE was first suggested by Tarasov and Zaslavsky [11,12] for fractal media. Since then, the fractional Ginzburg-Landau equation (FGLE) has been exploited to describe various physical phenomena, such as the media with fractal mass dimension [11], a class of critical phenomena when the organization of the system near the phase transition point is influenced by a competing nonlocal ordering [13] and a network of diffusively Hindmarsh-Rose neurons with long-range synaptic coupling [14]. Recently, the coupled fractional Ginzburg-Landau equation (CFGLE) with stocastic noise was discussed by [15].The dissipative mechanism of the FGLE and CFGLE are not characterized by the classical Laplacian but by the fractional power of the Laplacian, which brings some essential difficulties in theoretical analysis. Many authors have worked on the FGLE and CFGLE from the theoretical aspects [15][16][17][18][19][20]. For example, the solution of FGLE is derived by using psi-series with fractional powers [16]. Pu and Guo [17] investigated the global well-posedness and dynamics for the nonlinear FGLE. Millot and Sire [20] considered the asymptotic analysis of the FGLE in a bounded domain. Shu et al. [15] studied the random attractors for the CFGLE with stochastic noise. From the numerical point of view, there are quite a lot of numerical studies for the classical GLE, see [21-26] and references therein. But there are not too much numerical studies for the FGLE. Mvogo et al. [14] proposed a