# Optimized Computation of Tight Focusing of Short Pulses Using Mapping to Periodic Space
## 2 Materials and Methods
### Mapping to and from the Computational Subregion  
The stated problem of advancing electromagnetic field can be solved numerically by sampling and iterating the electromagnetic field in the region \((y\in(-R_{0},R_{0}))\cap(z\in(-R_{0},R_{0}))\), which has volume equal to approximately \(4R_{0}^{3}\) (assuming that \(R_{0}\gg L\)). However, geometric arguments indicate that the pulse should remain, approximately, in a thin spherically curved layer \(R_{0}-L/2-ct<R<R_{0}+L/2-ct\). Thus only a small fraction of the computational grid actually contains the pulse. If we were using the Finite Difference Time Domain (FDTD) method, we could reduce computational costs by updating the grid nodes in the nonempty regions only. However, as we wish to avoid numerical dispersion, we use a Fourier-based solver and thus do not have such a straightforward opportunity. We now describe how we can effectively exclude the empty regions and achieve a significant reduction of computational costs in this case. Let us consider a periodic field structure containing a series of repetitions of the pulse, spaced apart by a distance \(D\) along the \(x\) axis (see Figure 2a). If \(D\) is large enough, the pulses overlap neither initially nor during their propagation towards the focal region (see Figure 2b). When each pulse reaches its focus, a sequence of identical structures of focused pulses is formed (see Figure 2c). These structures are spaced apart by the distance \(D\) and have only a minor overlap due to diffractive spreading of the signal from the fringes of the initial structure. We can simulate this process by computing the evolution of the field in a region that has periodic boundary conditions and size \(D\) along \(x\) axis. The volume of this region is \(4 DR_{0}^{2}\ll 4R_{0}^{3}\) (we assume that \(D\sim L\ll R_{0}\)). Periodic boundary conditions are effectively provided by the Fourier solver by default. In this way we can obtain the numerical solution of the stated problem at reduced computational cost, with respect to both run time and memory. Let us describe particular details of our implementation. The choice of parameter \(D\) can be restricted by requiring that neighboring pulses do not overlap. From the geometrical consideration shown in Figure 3, one can obtain the following restriction:  
\[D\geq D_{min}=-R_{min}\cos\theta+\sqrt{R_{max}^{2}-R_{min}^{2}\sin^{2}\theta}, \tag{12}\]  
Figure 2: Illustration of the idea of the proposed method. The pictures show the evolution of the electromagnetic field composed by the equidistant replication of the initial pulse along the \(x\) axis: the field intensity is shown for (**a**) \(t=-R_{0}/c\), (**b**) \(t=-R_{0}/2c\), (**c**) \(t=0\) and (**d**) \(t=R_{0}/2c\) with f\(-\)number\(=0.3\) (\(\theta\approx 1\) rad), \(\lambda=1\) μm, \(R_{0}=16\lambda\), \(L=2\lambda\), \(\varepsilon=0.1\) rad, \(P_{0}=1\) W.  
where \(R_{max}=R_{0}+L/2\) and \(R_{min}=R_{0}-L/2\) are the boundaries of the initial pulse in spherical coordinates. Note that one can avoid this restriction by considering a superposition of the overlapping replicated pulses. In this case we can still obtain the separated structure of a single pulse at the focus if \(D>L\) due to the linearity of the Maxwell's equations. However, due to overlap, in this case we are not able to obtain (in a straightforward way) the intermediate states of the pulse on the way to the focus. The choice of \(D\) can also depend on the extent of fringes and the required accuracy. There is a trade-off: increasing \(D\) reduces the errors but increases the computational costs of the method. Although this choice of periodic domain is arbitrary with respect to translation along \(x\), it is useful to describe a particular implementation in detail. We suggest to carry out calculations in a layer \(x\in[-R_{min}-D,-R_{min}]\) with periodic boundary conditions.