  In this paper, we first develop a mathematical model for long-range,
hydrophobic attraction between amphiphilic particles. The non-pairwise
interactions follow from the first variation of a hydrophobic attraction domain
functional. The variation yields a hydrophobic stress that is used to
numerically calculate trajectories for a collection of two-dimensional
particles. The functional minimizer that accounts for hydrophobicity at
molecular-aqueous interfaces is a solution to a boundary value problem of the
screened Laplace equation. We reformulate the boundary value problem as a
second-kind integral equation (SKIE), discretize the SKIE using a Nystr\"om
discretization and `Quadrature by Expansion' (QBX) and solve the resulting
linear system iteratively using GMRES. We evaluate the required layer
potentials using the `GIGAQBX' fast algorithm, a variant of the Fast Multipole
Method (FMM), yielding the required particle interactions with asymptotically
optimal cost. The entire scheme is adaptive, high-order, and capable of
handling close-to-touching geometry. The simulated particle systems exhibit a
variety of multiscale behaviors over both time and length: Over short time
scales, the numerical results show self-assembly for model lipid particles. For
large system simulations, the particles form realistic configurations like
micelles and bilayers. Over long time scales, the bilayer shapes emerging from
the simulation appear to minimize a form of bending energy.
