I'm not a mathematical physicist, so parts of this may be wrong.
In quantum field theory, one encounters operators that are supported at points in spacetime, or at least are very local. For example (if we ignore uncertainty), one may have a "photon creation at $x$ with momentum $p$ and helicity $\xi$" operator, which changes the universe so that there is an extra photon at $x$ with momentum $p$ and helicity $\xi$. If you have encountered raising and lowering operators when studying the quantum mechanics of a harmonic oscillator, this is a generalization to the situation where you have harmonic oscillators at every point (and in interacting theories, the oscillators may coupled nonlinearly).
If we set up multiple operators supported at different points, there may be no canonical way to compose them. For example, if we have creations at spacelike separated points, any choice of time-ordering can be permuted by a boost. However, we can consider a parameter space of all possible ways to compose such operators, and we can consider how the compositions change when we continuously vary the positions of the creation events in spacetime. Furthermore, if we put a subset of such creations close to each other, we may view their composition as a single operator supported on the neighborhood, that depends on their relative position inside that neighborhood. Operads are precisely the objects that encode such families of composition laws and iterations.
One example where these arise is in the theory of vertex operator algebras, which (so I'm told) describe the chiral parts of conformal field theories. Given an compact complex algebraic curve and a vertex operator algebra, you can choose points with frames, and insert states (i.e., apply certain operators) at those points. There is a procedure by which you can construct from these data a space of chiral correlation functions, and this space will vary smoothly with the points and frames (and also with the complex moduli of the curve). If you put some points with frames close to each other and ignore the global geometry of the curve, the configurations of points with frames form a space in the framed $E_2$ operad, and the spaces of correlation functions are something like an algebra over the operad. This is a genus zero restriction of the compatibilities we see in conformal field theory, and can be strengthened to the statement (proved here under some extra hypotheses) that modules of a vertex operator algebra form a framed $E_2$ (in particular, braided) tensor category.