A mathematical formulation of the problem of electromagnetic wave scattering by a system of two homogeneous triaxial dielectric ellipsoids of complex index of refraction is presented. The analysis is based on the Lippman–Schwinger integral equation for an electric field. The corresponding integral equation for the scattering, which contains two singular kernels, is transformed into a pair of nonsingular integral equations for the angular Fourier transform of the electric field inside each scatterer. The latter equations are solved by reducing them by quadrature into a matrix equation. The resulting solutions are used to calculate the scattering amplitude. As a numerical application, the case of a two red blood cell rouleau model is considered. Typical values of the appropriate discretization parameters, which proved sufficient for achieving convergence, are presented, along with validity tests. The effect of the electromagnetic coupling of the scatterers is also illustrated. Efficient techniques, which are capable of reducing the rather high computing demands of the analysis, such as parallel processing, are both suggested and applied.