I'm asking about the solution of the Helmholtz equation on a periodic domain with piecewise constant wavespeed in different frequency regimes. One possible approach is to solving this problem is to write down integral equations on the boundary surfaces in terms of the Green's function of the system. Since the domain is periodic, this will be a periodic Green's function like $$ G(r,r') = \sum_L G_0(r,r'+L) $$ where $L$ is a lattice vector, and $G_0$ is something like $$ G_0(r,r') = \frac{e^{ik(r-r')}}{|r-r'|} $$ My question concerns how the computational cost of this method scales with frequency ($k$). Does it become computationally more difficult at low or high frequencies due to needing to include more terms in the lattice sum?
Edit: People seem to be answering different questions than the one I'm asking. I should clarify that I'm not interested in implementing such a method. I'm merely asking about the theoretical difficulties, as background for understanding the strengths and (moreso) weakness of the method in different frequency regimes. The kinds of problems I have in mind are (more or less) computing modes of periodic arrays of waveguides.