Many problems in Engineering and Biology necessitate solving the first passage time problem, which addresses questions such as the expected time for a Brownian particle in unbounded space to reach a target. I will present a boundary integral equation method for solving this mean first passage time with complex geometries of absorbing and reflecting bodies. The method applies the Laplace transform to the time-dependent problem, yielding a modified Helmholtz equation which is solved with a boundary integral method. This approach circumvents the limitations of traditional time-stepping methods and effectively handles the long equilibrium timescales associated with diffusion problems in unbounded domains. Returning to the time domain is achieved by applying quadrature along the so-called Talbot contour. I will demonstrate the method for various complex geometries formed by disjoint bodies of arbitrary shape on which either uniform Dirichlet or Neumann boundary conditions are applied. The examples include geometries that guide diffusion processes to particular absorbing sites, absorbing sites that are shielded by reflecting bodies, and finding the exits of confining geometries, such as mazes.