<p>This work investigates a three-dimensional coupled system of distributed-order fractional differential equations involving Caputo-type time derivatives and Riesz fractional spatial operators. The model is designed to capture complex physical phenomena characterized by memory effects and spatial nonlocality, which commonly arise in applications such as porous media flow, viscoelastic materials, and anomalous diffusion processes. To efficiently solve this system, we develop a fully discrete numerical scheme that integrates a high-order temporal discretization based on quadrature approximations of the distributed-order Caputo derivatives with a meshless collocation method for spatial discretization. The proposed approach offers flexibility in handling irregular geometries without the need for structured grids. Rigorous convergence and stability analyses are presented, yielding optimal error estimates under suitable regularity assumptions. Numerical experiments are performed to validate the accuracy and efficiency of the method, including tests with both exact solutions and more realistic scenarios lacking analytical solutions. The results demonstrate that the proposed scheme is accurate, stable, and computationally efficient, making it well suited for practical applications involving high-dimensional, nonlocal fractional systems.</p>