## Abstract

Wavefront shaping is a powerful technique that can be used to focus light through scattering media, which can be important for imaging through scattering samples such as tissue. This method is based on the assumption that the field at the output of the medium is a linear superposition of the modes traveling through different paths in the medium. However, when the scattering medium also exhibits nonlinearity, as may occur in multiphoton microscopy, this assumption is violated and the applicability of wavefront shaping becomes unclear. Here, using a simple model system with a scattering layer followed by a nonlinear layer, we show that with adaptive optimization of the wavefront, light can still be controlled and focused through a scattering medium in the presence of nonlinearity. Notably, we find that moderate positive nonlinearity can serve to significantly increase the focused fraction of power, whereas negative nonlinearity reduces it.

© 2017 Optical Society of America

## 1. INTRODUCTION

A major limitation of optical imaging is the inability to image deep into inhomogeneous media. Inhomogeneity causes scattering, which randomizes the direction of propagation of the light and prevents focusing or imaging with a lens. Recently, Vellekoop and Mosk have demonstrated that optimization of the incident wavefront using a spatial light modulator (SLM) can be used to partially restore the diffraction-limited focus [1]. Conceptually, each pixel of the SLM can be viewed as the source of a different speckle pattern, and by devising a clever phase mask, all these patterns can be made to constructively interfere at a given point at the output plane, enhancing its intensity by a factor that is approximately equal to the number of pixels. This successful demonstration prompted a considerable amount of research on the feasibility and applicability of wavefront shaping in many fields [2–10], including imaging through scattering tissue [11–13]. One challenge in the application of this technique to imaging through scattering media is that, while the brightness of the focus can be increased by 2–3 orders of magnitude, the power contained in the obtained focus is only a small fraction of the total available power. This limitation is imposed by the finite number of degrees of control that can be utilized in practice in wavefront shaping, which is typically significantly smaller than the number of transmission modes through the scattering medium, even when it is thin [14]. The rest of the light remains in the form of a speckle pattern, which creates a strong background signal and may completely overwhelm the signal from the obtained focus. Additionally, since the power at the focus is limited, the total intensity of the incident light must be significantly increased to compensate and may damage the sample.

When the medium exhibits nonlinearity in addition to inhomogeneity, the field at the output plane can no longer be described as a linear superposition of speckle fields originating from the SLM pixels, since these fields interact as they propagate [15,16]. As a result, commonly used methods of focusing through scattering media may become ineffective. Transmission-matrix-based methods, for example, rely on the linear superposition principle and therefore may be inadequate. Furthermore, time reversal of the scattered wave is highly challenging since any inaccuracy in the reconstructed beam is amplified by nonlinear propagation. Determining if and to what degree the propagation of light can still be controlled in scattering media with nonlinearity is both a fundamental and practical question, with implications for techniques employing short pulses of light, including nonlinear imaging, laser microsurgery, and nonlinear photodynamic therapy [17–22]; focusing of light through multimode fibers [23–25]; as well as for applications employing high-power beams, which experience thermal nonlinearity in a variety of samples [26,27].

In this work, we study the focusing of light through a scattering medium in the presence of nonlinearity using adaptive optimization of the input wavefront, and show that focusing nonlinearity in fact serves to increase the fraction of total power that the focus contains. Since adaptive optimization iteratively adjusts the solution as opposed to relying on calibration measurements, it is well suited for nonlinear media. To elucidate the roles of different elements in the scheme, we studied wavefront shaping in a simple model system consisting of a thin scattering layer followed by a non-scattering nonlinear layer, both numerically and experimentally. We show numerically that this simplified system produces the same qualitative results as a system in which the scattering and nonlinearity are distributed throughout the medium.

We note that speckles propagating through nonlinear media often change their width and shape due to self-focusing or defocusing effects [28]. This change affects the maximum intensity of the speckle without affecting the power it contains. As a result, examining the enhancement of a given speckle by wavefront shaping, defined as the ratio of the maximum intensity of the enhanced speckle to the mean intensity of the unoptimized speckle field, can be ambiguous in nonlinear media. Therefore, throughout this work, we will analyze the *focused power fraction* achieved by wavefront shaping, which is the ratio of the power contained in the enhanced speckle, or focus, to the total power contained in the speckle field.

Naively, one may expect that the combination of wavefront shaping and nonlinearity can be thought of as two separate effects operating on the speckle pattern sequentially. That is, the total effect can be predicted by combining the increase in power in the focus expected due to wavefront shaping in a linear medium, with the expected self-focusing or self-defocusing of a single speckle caused by the nonlinear medium. Consequently, the focus will become narrower (broader) in the presence of focusing (defocusing) nonlinearity due to the self-action of the enhanced speckle, but the fraction of incident power contained in this speckle will not change compared to the linear case. In fact, our results indicate that the combination of the two effects leads to a modified collective behavior, which, for focusing nonlinearity, allows for a greater fraction of the incident power to be controlled and focused. The opposite is true for defocusing nonlinearity, where the fraction of incident power that can be controlled is smaller than for a linear medium. To explain this result, we will first discuss how nonlinearity modifies the properties of speckle fields in general before we consider the added effect of wavefront shaping.

## 2. NONLINEAR PROPAGATION OF 2D SPECKLE FIELDS

The propagation of a two-dimensional speckle field inside a medium with Kerr-type nonlinearity is described by the nonlinear wave equation [29]

Notably, we found that the change in speckle statistics is caused in part by a redistribution of the total power between the speckles and not just by the narrowing (broadening) of the speckles themselves, which would result only in higher (lower) peak intensities. Figure 1(e) depicts the variation of the *total* power contained in the brightest speckle in the speckle field after propagation through media with varying amounts of nonlinearity, averaged over many realizations (*without* wavefront shaping). The nonlinearity values are given as a normalized quantity, $\eta \equiv {I}_{\text{in}}{k}_{0}{n}_{2}L$, where ${I}_{\text{in}}$ is the intensity of the incident field and $L$ is the total propagation distance. $\eta $ therefore represents the total nonlinear phase accumulated during propagation in the medium. The two extreme points in Fig. 1(e), namely $\eta =-0.66$ and $\eta =0.42$, match the values of nonlinearity used to create Figs. 1(b) and 1(c), respectively. The power contained in the speckle was estimated by fitting it to a 2D Gaussian function and integrating the result. The values of power are given relative to the power contained in the brightest speckle in a linear medium ($\eta =0$). We see that the power contained in the speckle declines with $\eta $ for defocusing nonlinearity and grows with $\eta $ for focusing nonlinearity. Hence, focusing nonlinearity redistributes the total power among the speckles such that a smaller number of bright speckles holds a larger fraction of the total power of the field. Defocusing nonlinearity redistributes the total power such that it is spread out more evenly, among many speckles. We expect that this redistribution of power may influence the fraction of incident power that can be controlled and focused using wavefront shaping.

## 3. SIMULATION RESULTS

Our prediction was tested by performing simulations of wavefront shaping of the field incident upon a forward-scattering layer (diffuser) followed by a nonlinear layer. Between the scattering and nonlinear layers, the field was allowed free-space propagation so that the speckle field entering the nonlinear layer was fully developed. The spatial phase of the field incident upon the diffuser was optimized adaptively in order to enhance a single speckle at the output of the nonlinear medium (the configuration of the simulation follows the experimental setup, depicted in Fig. 4). This optimization was repeated for different input powers, corresponding to different nonlinearity strengths. The results of the simulation are shown in Fig. 2. The data was fitted (black lines) by simple functions of exponential form as a guide to the eye (see Section 6, Methods). We can see that the fraction of power that can be controlled and focused indeed decreases [blue circles in Fig. 2(a)] or increases [red circles in Fig. 2(b)] significantly as a function of nonlinearity strength. This result shows that wavefront shaping causes a dynamic redistribution of power between the speckles. In particular, in a focusing nonlinear medium, there is a positive feedback effect: as power is being focused into the region of the chosen speckle by wavefront shaping, $\eta $ increases and the power becomes redistributed such that a smaller amount of speckle holds a larger fraction of the total power. This allows for more power to be controlled and focused by wavefront shaping, which increases $\eta $, and so on. The positive feedback leads to a focus with a larger fraction of the total power, even for $\eta $ values that are relatively low for the unshaped field. For example, from Fig. 2(b) we can see that the power obtained in the enhanced speckle (the focus) is twice as large for $\eta =0.062$ than for the linear case. Yet from Fig. 1(e) we see that, for $\eta \simeq 0.06$, the distribution of power between the speckles in the initial unshaped field is approximately the same as for the linear case. In a defocusing nonlinear medium, however, as more power is focused into a certain region, it spreads out more evenly between a large number of speckles, hindering the ability to focus power to a chosen speckle at the output. The result is a focus that contains a smaller fraction of the total power than for a linear medium.

Moreover, we propose that our simplified model, consisting of a scattering layer followed by a nonlinear layer, exhibits the same qualitative behavior as a scattering nonlinear medium, in which scattering and nonlinearity are homogeneously distributed. In order to test this assumption, we performed simulations of focusing coherent light through a composite medium, consisting of multiple alternating scattering and nonlinear layers [34]. The configuration used for the multi-layer simulation is identical to that of the single-layer simulation except for the modification of the medium, and that there is no free-space propagation between the scattering and nonlinear layers. A sketch of the configuration used is shown in Fig. S1 of Supplement 1. The results of the composite medium simulation are presented in Fig. 2(c). The focused power fraction as a function of normalized nonlinearity $\eta \equiv {I}_{\text{in}}{k}_{0}{n}_{2}L$ is represented by red circles. $L$ in this case is the total length of all nonlinear layers. The focused power fractions are normalized to that obtained in a linear multi-layer medium ($\eta =0$). The results follow the same trend as those of the single-layer simulation [Fig. 2(a)], showing focused power fractions that grow by approximately a factor of 2 as focusing nonlinearity is increased (the $\eta $ values are larger for this simulation due to the removal of the free-space propagation; see Supplement 1). Therefore, we may conclude that the trends we observe with our simplified model, consisting of a single scattering layer followed by a single nonlinear layer, are general and hold also for a system with scattering and nonlinearity distributed throughout the medium.

The mechanism responsible for the increase in the focused power can be explored by examining the propagation of the optimized field through the nonlinear layer, in the single-layer simulation. Figure 3 shows a one-dimensional cut through the simulated two-dimensional speckle field after wavefront optimization, as it propagates from the input facet of the nonlinear layer (propagation $\text{length}=0$) towards its output (propagation $\text{length}=1$), where a focus is created. Since the three fields were optimized to focus after propagation in layers with different nonlinearity types, the optimal SLM phase is different in all three cases and thus also the speckle field entering the nonlinear layer, at propagation $\text{length}=0$ (see Fig. S2 in Supplement 1 for a version of Fig. 3 without wavefront shaping). In both nonlinear cases, the propagation of the speckle field is dynamic, with bright speckles changing quickly into dark and vice versa, whereas in the linear case the propagation is quite static. Accordingly, in the linear case, the wavefront optimized to create a focus at the output already contains a bright speckle in that location at the input, which develops to half its final intensity halfway through the medium [Fig. 3(a)]. In the nonlinear cases, both focusing [Fig. 3(b)] and defocusing [Fig. 3(c)], the optimized wavefront does not initially contain a bright speckle at the focus location. Furthermore, the enhanced speckle reaches half of its maximal intensity only close to the output, at propagation length $\simeq 0.9$, creating a tightly confined focus along the propagation direction. In addition, one can observe in Fig. 3(c) that several bright speckles follow trajectories leading into the enhanced speckle as they propagate from the input to the output, almost as we would expect when focusing a non-speckled field with a lens (or with self-focusing). Thus, in a focusing nonlinear medium, wavefront control can be used to direct energy from neighboring bright speckles into the enhanced speckle.

## 4. EXPERIMENTAL RESULTS

#### A. Focusing Light Through Scattering Media in the Presence of Nonlinearity

Our experimental scheme is presented in Fig. 4. Coherent light was passed through an SLM, where it acquired a spatially dependent phase. The shaped light was imaged with a telescope onto a thin scattering medium and, after some free-space propagation, entered a nonlinear medium where it propagated through 5 cm. The output facet of the nonlinear medium was imaged onto a CCD camera. An iterative algorithm searched for the optimal SLM phase, which maximized the power of a chosen speckle at the output of the nonlinear medium. For the focusing nonlinearity experiment, the light source was a pulsed femtosecond laser and the nonlinear medium used was ethanol, which exhibits Kerr effect. For the defocusing nonlinearity experiment, the light source was a CW laser and the nonlinear medium used was a weakly absorbing dye diluted in ethanol, which exhibits defocusing (thermal) nonlinearity [28,30] (for more details, see Section 6, Methods).

The experimental results of the obtained focused power fractions, at different laser powers, are presented in Fig. 5. Figure 5(a) shows the fraction of focused power in a defocusing nonlinear medium (blue circles) as a function of the input laser power, and the scaled simulation results for comparison (black line). Figure 5(b) shows the fraction of focused power in a focusing nonlinear medium (red circles) as a function of the peak power incident on the diffuser, as well as the scaled simulation results for comparison (black line). The focused power fractions in both figures are relative to the linear case. In both experiments, the results follow the trend of the simulation, either increasing or decreasing with nonlinearity strength. The experimental results thus verify the prediction that moderate nonlinearity significantly alters the fraction of power that can be controlled and focused through a scattering medium, and that the focused power fraction can be substantially increased in the presence of mild focusing nonlinearity. The range of powers shown in Fig. 5(a) corresponds to the low nonlinearity range of the simulation [$\eta \lesssim 0.28$ in Fig. 2(a)], while the range of powers shown in Fig. 5(b) approximately corresponds to the full nonlinearity range of the simulation [Fig. 2(b)]. The peak powers used in the focusing nonlinear experiment are over an order of magnitude lower than those typically used at the focus during two-photon microscopy, which are approximately in the $\sim 10\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$ range. Therefore, even taking into account the $\sim 20$-fold increase in $\eta $ values predicted by the multi-layer simulation, we expect significant changes in the focused power fraction to occur in scattering samples with high water content during nonlinear microscopy measurements. We note that the nonlinearity of the defocusing medium used is fundamentally nonlocal, whereas the simulated defocusing nonlinearity is local. Yet the agreement between the experimental and simulation results shows that the nonlocality does not significantly alter the results of our experiment.

#### B. “Nonlinear Memory Effect”

Another interesting aspect of focusing of light through nonlinear scattering media is the robustness of the optimization to fluctuations or changes in the nonlinearity strength. Such fluctuations can be caused, for example, by fluctuations in the intensity of the laser. To address this question experimentally, after finding the optimal SLM phase for a certain laser power, we measured the fraction of focused power obtained with this SLM phase for other laser powers. This measurement was performed for all optimization powers. In a sense, this is a measure of the “nonlinear memory effect” [35] of the scattering medium, i.e., the degree to which the optimal SLM phase can be used to obtain a focus as the nonlinearity strength of the medium is varied, without reoptimization.

The results are presented in Fig. 5(c) for defocusing nonlinearity and Fig. 5(d) for focusing nonlinearity. The different colors represent measurements performed for different optimization powers, shown in the legend. The results for the two cases are qualitatively different. In the defocusing case, we see a clear drop in the focused power on either side of the optimization power, and the drop is roughly symmetrical between the two sides. Furthermore, for larger optimization powers corresponding to larger nonlinearity values (and lower focused power fractions), the widths of the curves grow. Therefore, the optimization is more robust and the optimal SLM phase can be used to obtain a focus of similar quality for larger variations in the nonlinearity strength compared to lower nonlinearity values. In contrast, in the focusing nonlinear case, the curves for all optimization powers are quite similar (the curves of the other four optimization powers were omitted from the figure for clarity). The optimization is very robust and can be used to obtain a focus of similar quality all throughout the nonlinearity range examined in the experiment ($0.1\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}<\text{laser}$ peak power $<0.5\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$) and beyond it, up to $\sim 1\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$. For powers below $0.1\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$, the focused power fraction drops rapidly, whereas for power above $1\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$ the focused power fraction declines more slowly, reaching half of the original value only around $\sim 5\u201310$ times the optimization powers.

## 5. DISCUSSION

We have demonstrated, both numerically and experimentally, focusing of light through a scattering layer followed by a nonlinear layer, and shown that nonlinearity significantly modifies the fraction of power that can be channeled to the focus. In particular, the presence of moderate focusing nonlinearity has been shown to cause about a two-fold increase in the total power in the obtained focus compared to a linear medium. Conversely, the presence of moderate defocusing nonlinearity has been shown to cause approximately a two-fold decrease in the focused power fraction. Our simulations show that these findings are valid also when the scattering and nonlinearity are homogeneously distributed throughout the medium.

The results presented in this work suggest the favorability of focusing and imaging through scatterers with short pulses of light, which experience mild Kerr self-focusing in a variety of samples, including biologically relevant samples [36–38]. Previous works have shown that pulsed fields that propagate through scattering layers maintain temporal coherence long after the spatial coherence is lost [7], and that nonlinear signals can be utilized to focus light noninvasively [39]. Our current results show that the basic premises of wavefront shaping indeed extend into the nonlinear regime, and we expect that such noninvasive nonlinear techniques might even be enhanced by the nonlinear response of the specimen.

## 6. METHODS

#### A. Single-Layer Simulation Implementation

Nonlinear propagation was simulated using the split-step method. Time-domain was not addressed specifically in the simulation, as this would result in impractical run times. The diffuser was modeled as a single layer, with no thickness, of randomly distributed phase features $\varphi $ in the range of $-\pi <\varphi <\pi $. The scattering is therefore only in the forward direction.

Defocusing nonlinearity was modeled as Kerr-type, with the nonlinear operator $\widehat{N}=-{k}_{0}{n}_{2}{|E|}^{2}$. The beam illuminated $\sim 110\times 110$ features on the modeled diffuser, which resulted in a speckle field that contained several thousands of speckles. The SLM was modeled with $15\times 15\text{\hspace{0.17em}}\mathrm{pixels}$. The linear ($\eta =0$) enhancement value and focused power fraction were 300 and 0.022, respectively. However, in the simulation *without* wavefront shaping shown in Fig. 1(b), the implementation of the defocusing simulation matched that of the focusing simulation exactly, for the purpose of comparison. Therefore, the beam illuminated $\sim 30\times 30$ features on the modeled diffuser, which resulted in a speckle field that contained several hundreds of speckles.

Focusing nonlinearity was modeled with a higher-order defocusing term, $\widehat{N}={k}_{0}({n}_{2}{|E|}^{2}-{n}_{4}{|E|}^{4})$, where ${n}_{4}\ll {n}_{2}$ for the simulations with and without wavefront shaping. The higher-order term introduces saturation of the self-focusing process, as is often observed in various Kerr samples [40], and prevents the collapse of the speckles to sizes below the resolution of the simulation. The beam illuminated $\sim 30\times 30$ features on the modeled diffuser, which resulted in a speckle field that contained several hundreds of speckles. The SLM was modeled with $8\times 8\text{\hspace{0.17em}}\mathrm{pixels}$. The simulation parameters were modified for the focusing case because of its increased run time (due to stronger nonlinearity effects that require more propagation steps). The simulated window size was therefore reduced in order to compensate for the increase. In order to equate the nonlinearity scales of the two simulations, the $\eta $ values for the defocusing simulation were divided by a factor of 2. The linear ($\eta =0$) enhancement value and focused power fraction were 77 and 0.0413, respectively.

For both simulations, the speckles entering the media were fully developed, and $\eta $ was varied by varying ${I}_{\text{in}}$. An arbitrary speckle at the output of the nonlinear layer was chosen and optimized using an iterative genetic algorithm. The functions used for fitting the simulation results in Fig. 2 were: (a) $f(x)=\frac{a}{\mathrm{exp}(bx-c)+1}+d$ and (b) $f(x)=a\text{\hspace{0.17em}}\mathrm{exp}(-bx)+c$. For details on the implementation of the multi-layer simulation, see Supplement 1.

#### B. Experimental Implementation

In the defocusing nonlinearity experiment, light from a 532 nm CW laser was expanded with a $6\times $ magnifying telescope and reflected off a two-dimensional liquid-crystal SLM (Hamamatsu LCOS X10468-01). The phase-shaped wavefront was imaged with a $6\times $ demagnifying telescope onto a holographic diffuser (5° diffusion FWHM) and propagated $\sim 5\text{\hspace{0.17em}}\mathrm{cm}$ in air before it entered the nonlinear medium. The transmission of the diffuser was $\sim 90\%$. The output facet of the nonlinear medium was then imaged onto a CCD camera (Andor Luca S). An arbitrary speckle was chosen and optimized using an iterative genetic algorithm. The nonlinear medium was LDS 751 dye dissolved in ethanol contained in a 5 cm path-length cylindrical cuvette. The solution weakly absorbed the light due to the low concentration of the dye, and exhibited thermal nonlinearity due to the ethanol. Of the total laser power, $\sim 70\%$ entered the cuvette and $\sim 50\%$ of that was absorbed in it. The genetic algorithm was terminated for all measurements after 1550 generations, a number of generations that showed sufficient convergence for all runs. The enhancement value for the linear run (0.5 mW input power) was 650. An accurate absolute measure of the focused power fraction is more challenging to obtain in the experimental implementation due to the finite aperture of the CCD, and therefore only relative focused power fractions were calculated.

In the focusing nonlinearity experiment, 35 fs pulses from an amplified femtosecond laser (bandwidth of $\sim 30\text{\hspace{0.17em}}\mathrm{nm}$, centered at 800 nm) were significantly attenuated using a waveplate and a polarizing beam splitter to obtain the desired peak powers. The beam was directly reflected off a broadband two-dimensional SLM (Hamamatsu LCOS X10468-02). The phase-shaped wavefront was imaged with a $3.33\times $ demagnifying telescope onto a holographic diffuser (1° diffusion FWHM) and propagated $\sim 15\text{\hspace{0.17em}}\mathrm{cm}$ in air before it entered the nonlinear medium. The output facet of the nonlinear medium was then imaged onto the CCD camera, where a chosen speckle was optimized using an iterative genetic algorithm. The nonlinear medium was neat ethanol contained in a 5 cm path-length cylindrical cuvette. The peak powers provided in Figs. 5(b) and 5(d) are those incident on the diffuser. The genetic algorithm was terminated for all measurements after $\sim 1300$ generations, which showed sufficient convergence for all runs. The enhancement value for the linear run ($3.6\mathrm{e}\text{-}3\text{\hspace{0.17em}}\mathrm{GW}/{\mathrm{cm}}^{2}$ incident peak power) was 300.

#### C. Evaluation of Focused Power

The amount of power contained in the focused speckle, in both the simulations and the defocusing nonlinearity experiment, was calculated by fitting it with a two-dimensional Gaussian function with a different width in each dimension and integrating the fitted expression. The Gaussian fit was not accurate enough in the focusing nonlinearity experiment, and therefore the energy in the focused speckle was calculated simply by integration over the speckle area. The error bars for the $\eta =0$ data point for Figs. 2 and 5 are presented as an example of the typical error in these plots. The error values were calculated by applying the focusing procedure 20 times using the same parameters and calculating the standard deviation of the obtained focused power fractions for the simulations and defocusing experiment. For the focusing experiment, the errors were calculated with 9 runs.

## Funding

Israel Science Foundation (ISF) (Icore on Light and Matter).

## Acknowledgment

The authors thank R. Fischer, M. Segev, H. H. Sheinfux, and D. Gilboa for helpful discussions. The authors gratefully acknowledge funding from the ISF.

See Supplement 1 for supporting content.

## REFERENCES

**1. **I. M. Vellekoop and A. P. Mosk, “Focusing coherent light through opaque strongly scattering media,” Opt. Lett. **32**, 2309–2311 (2007). [CrossRef]

**2. **T. Cizmar, M. Mazilu, and K. Dholakia, “In situ wavefront correction and its application to micromanipulation,” Nat. Photonics **4**, 388–394 (2010). [CrossRef]

**3. **I. M. Vellekoop, A. Lagendijk, and A. P. Mosk, “Exploiting disorder for perfect focusing,” Nat. Photonics **4**, 320–322 (2010). [CrossRef]

**4. **S. Popoff, G. Lerosey, M. Fink, A. C. Boccara, and S. Gigan, “Image transmission through an opaque material,” Nat. Commun. **1**, 81 (2010). [CrossRef]

**5. **E. G. van Putten, D. Akbulut, J. Bertolotti, W. L. Vos, A. Lagendijk, and A. P. Mosk, “Scattering lens resolves sub-100 nm structures with visible light,” Phys. Rev. Lett. **106**, 193905 (2011). [CrossRef]

**6. **X. Xu, H. Liu, and L. V. Wang, “Time-reversed ultrasonically encoded optical focusing into scattering media,” Nat. Photonics **5**, 154–157 (2011). [CrossRef]

**7. **O. Katz, E. Small, Y. Bromberg, and Y. Silberberg, “Focusing and compression of ultrashort pulses through scattering media,” Nat. Photonics **5**, 372–377 (2011). [CrossRef]

**8. **J. Aulbach, B. Gjonaj, P. M. Johnson, A. P. Mosk, and A. Lagendijk, “Control of light transmission through opaque scattering media in space and time,” Phys. Rev. Lett. **106**, 103901 (2011). [CrossRef]

**9. **D. J. McCabe, A. Tajalli, D. R. Austin, P. Bondareff, I. A. Walmsley, S. Gigan, and B. Chatel, “Spatio-temporal focusing of an ultrafast pulse through a multiply scattering medium,” Nat. Commun. **2**, 447 (2011). [CrossRef]

**10. **Y. Choi, T. D. Yang, C. Fang-Yen, P. Kang, K. J. Lee, R. R. Dasari, M. S. Feld, and W. Choi, “Overcoming the diffraction limit using multiple light scattering in a highly disordered medium,” Phys. Rev. Lett. **107**, 023902 (2011). [CrossRef]

**11. **Z. Yaqoob, D. Psaltis, M. S. Feld, and C. Yang, “Optical phase conjugation for turbidity suppression in biological samples,” Nat. Photonics **2**, 110–115 (2008). [CrossRef]

**12. **Y. M. Wang, B. Judkewitz, C. A. Dimarzio, and C. Yang, “Deep-tissue focal fluorescence imaging with digitally time-reversed ultrasound-encoded light,” Nat. Commun. **3**, 928 (2012). [CrossRef]

**13. **O. Katz, E. Small, and Y. Silberberg, “Looking around corners and through thin turbid layers in real time with scattered incoherent light,” Nat. Photonics **6**, 549–553 (2012). [CrossRef]

**14. **A. P. Mosk, A. Lagendijk, G. Lerosey, and M. Fink, “Controlling waves in space and time for imaging and focusing in complex media,” Nat. Photonics **6**, 283–292 (2012). [CrossRef]

**15. **C. Barsi, W. Wan, and J. W. Fleischer, “Imaging through nonlinear media using digital holography,” Nat. Photonics **3**, 211–215 (2009). [CrossRef]

**16. **A. Goy and D. Psaltis, “Digital reverse propagation in focusing Kerr media,” Phys. Rev. A **83**, 031802 (2011). [CrossRef]

**17. **D. Yelin and Y. Silberberg, “Laser scanning third-harmonic-generation microscopy in biology,” Opt. Express **5**, 169–175 (1999). [CrossRef]

**18. **C. L. Evans and X. Sunney Xie, “Coherent anti-Stokes Raman scattering microscopy: chemical imaging for biology and medicine,” Annu. Rev. Anal. Chem. **1**, 883–909 (2008). [CrossRef]

**19. **B. G. Saar, C. W. Freudiger, J. Reichman, C. M. Stanley, G. R. Holtom, and X. Sunney Xie, “Video-rate molecular imaging in vivo with stimulated Raman scattering,” Science **330**, 1368–1370 (2010). [CrossRef]

**20. **C. L. Hoy, O. Ferhanoglu, M. Yildirim, K. H. Kim, S. S. Karajanagi, K. M. C. Chan, J. B. Kobler, S. M. Zeitels, and A. Ben-Yakar, “Clinical ultrafast laser surgery: recent advances and future directions,” IEEE J. Quantum Electron. **20**, 7100814 (2014). [CrossRef]

**21. **A. Heisterkamp, T. Ripken, T. Mamom, W. Drommer, H. Welling, W. Ertmer, and H. Lubatschowski, “Nonlinear side effects of fs pulses inside corneal tissue during photodisruption,” Appl. Phys. B **74**, 419–425 (2002). [CrossRef]

**22. **J. D. Bhawalker, N. D. Kumar, C.-F. Zhao, and P. N. Prasad, “Two-photon photodynamic therapy,” J. Clin. Laser Med. Surg. **15**, 201–204 (1997).

**23. **S. Bianchi and R. Di Leonardo, “A multi-mode fiber probe for holographic micromanipulation and microscopy,” Lab Chip **12**, 635–639 (2012). [CrossRef]

**24. **T. Cizmar and K. Dholakia, “Shaping the light transmission through a multimode optical fibre: complex transformation analysis and applications in biophotonics,” Opt. Express **19**, 18871–18884 (2011). [CrossRef]

**25. **L. G. Wright, D. N. Christodoulides, and F. W. Wise, “Controllable spatiotemporal nonlinear effects in multimode fibres,” Nat. Photonics **9**, 306–310 (2015). [CrossRef]

**26. **M. Motamedi, A. J. Welch, W.-F. Cheong, S. A. Ghaffari, and O. T. Tan, “Thermal lensing in biologic medium,” IEEE J. Quantum Electron. **24**, 693–696 (1988). [CrossRef]

**27. **W.-C. Lin, M. Motamedi, and A. J. Welch, “Dynamics of tissue optics during laser heating of turbid media,” Appl. Opt. **35**, 3413–3420 (1996). [CrossRef]

**28. **Y. Bromberg, Y. Lahini, E. Small, and Y. Silberberg, “Hanbury Brown and Twiss interferometry with interacting photons,” Nat. Photonics **4**, 721–726 (2010). [CrossRef]

**29. **R. W. Boyd, *Nonlinear Optics*, 3rd ed. (Academic, 2008).

**30. **S. Derevyanko and E. Small, “Nonlinear propagation of an optical speckle field,” Phys. Rev. A **053816**, 053816 (2012). [CrossRef]

**31. **C. Sun, S. Jia, C. Barsi, S. Rica, A. Picozzi, and J. W. Fleischer, “Observation of the kinetic condensation of classical waves,” Nat. Phys. **8**, 471–475 (2012). [CrossRef]

**32. **R. E. de Araujo and A. S. L. Gomes, “Nonlinear optical Kerr coefficients of disordered media,” Phys. Rev. A **57**, 2037–2040 (1998). [CrossRef]

**33. **S. E. Skipetrov and R. Maynard, “Instabilities of waves in nonlinear disordered media,” Phys. Rev. Lett. **85**, 736–739 (2000). [CrossRef]

**34. **Y. Wang and J. W. Fleischer, “A model of thick scattering media using a series of random diffraction gratings,” in *Imaging and Applied Optics* (2017), paper JTu5A.

**35. **I. Freund, M. Rosenbluh, and S. Feng, “Memory effects in propagation of optical waves through disordered media,” Phys. Rev. Lett. **61**, 2328–2331 (1988). [CrossRef]

**36. **B. A. Rockwell, W. P. Roach, M. E. Rogers, M. W. Mayo, C. A. Toth, C. P. Cain, and G. D. Noojin, “Nonlinear refraction in vitreous humor,” Opt. Lett. **18**, 1792–1794 (1993). [CrossRef]

**37. **E. O. Potma, W. P. De Boeij, and D. A. Wiersma, “Femtosecond dynamics of intracellular water probed with nonlinear optical Kerr effect microspectroscopy,” Biophys. J. **80**, 3019–3024 (2001). [CrossRef]

**38. **N. T. Hunt, L. Kattner, R. P. Shanks, and K. Wynne, “The dynamics of water-protein interaction studied by ultrafast optical Kerr-effect spectroscopy,” J. Am. Chem. Soc. **129**, 3168–3172 (2007). [CrossRef]

**39. **O. Katz, E. Small, Y. Guan, and Y. Silberberg, “Noninvasive nonlinear focusing and imaging through strongly scattering turbid layers,” Optica **1**, 170–174 (2014). [CrossRef]

**40. **C. Bree, A. Demircan, and G. Steinmeyer, “Saturation of the all-optical Kerr effect,” Phys. Rev. Lett. **106**, 183902 (2011). [CrossRef]