The deformation and stress fields due to a three-dimensional, pressurized magma chamber are computed using the Indirect Boundary Integral Method (IBIM) with a numerical scheme based on point, single-force distribution over the closed surface of the chamber, and Green's function representation of the contribution of each single-force to the overall deformation. This scheme follows on Yang et al. (1988) analytical solution of the deformation due to ellipsoidal cavities, obtained from prescribed distributions of dilatation dipoles and double-couples in and around the chamber, respectively. Our solutions for generic cases of three-dimensional magma chambers in a half-space provide excellent agreement with those from analytical or other numerical, well accepted solutions. As well, our algorithm computes the deformation derivatives and stress fields via Hooke's law. Using the stress tensor field, we compute the Coulomb Failure Stress (CFS) on a target earthquake fault plane. We have calculated the CFS values caused by inflating spherical, prolate ellipsoidal, penny shaped and dyke chambers on a nearby fault. For all these cases, we prescribe the elastic parameters, the fault dimensions, orientation and dip, and the volumetric change of the magma chamber upon inflation. A parametric study in terms of the magma chamber dimensions shows that the CFS depends strongly on the shape and orientation of the chamber with respect to the fault. Although we do not present the results for arbitrarily shaped magma chambers other than those described above, we emphasize that the numerical formulation of our algorithm is not constrained by a particular shape, and can be applied to realistic cases of magma reservoirs. The only assumption is that the distributed single forces are perpendicular to the internal wall of the chamber, regardless of its shape. The magnitude of the forces is determined by the elastic properties of the confining rock, upon volumetric changes. The method is remarkably fast, judging from our solution of the benchmark case of an ellipsoidal magma chamber in a half-space. Using 1150 point forces to represent the chamber, it takes 5.54s CPU time to compute the ground deformation at 2500 points on the free-surface. This, we believe, can be useful in applications where the computation of many trial models are needed to converge to a final solution. Such applications could be the interpretation of surface deformation from GNSS data, and the study of volcano-fault interactions by computing the stress perturbations of volcanic episodes and nearby earthquake fault ruptures. [ABSTRACT FROM AUTHOR]