This manuscript introduces a discrete technique to estimate the solution of a doublefractional two-component Bose–Einstein condensate. The system consists of two coupled nonlinear
parabolic partial differential equations whose solutions are two complex functions, and the spatial
fractional derivatives are interpreted in the Riesz sense. Initial and homogeneous Dirichlet boundary
data are imposed on a multidimensional spatial domain. To approximate the solutions, we employ a
finite difference methodology. We rigorously establish the existence of numerical solutions along with
the main numerical properties. Concretely, we show that the scheme is consistent in both space and
time as well as stable and convergent. Numerical simulations in the one-dimensional scenario are
presented in order to show the performance of the scheme. For the sake of convenience, A MATLAB
code of the numerical model is provided in the appendix at the end of this work.
