Open Access
Add to BibSonomy
Abstract
The investigations on coherent enhancement absorption (CEA) inside scattering media are critically important in biophotonics. CEA can deliver light to the targeted position, thus enabling deep-tissue optical imaging by improving signal strength and imaging resolution. In this work, we develop a numerical framework that employs the method of finite-difference time-domain. Both the transmission and reflection matrices of scattering media with open boundaries are constructed, allowing the studies on the eigenvalues and eigenchannels. To realize CEA for scattering media with local absorption, we develop a genetic-algorithm-assisted numerical model. By minimizing the total transmittance and reflectance simultaneously, different realizations of CEA are observed and, without setting internal monitors, can be differentiated with cases of light leaked from sides. By modulating the incident wavefront at only one side of the scattering medium, it is shown that for a 5-μm-diameter absorber buried inside a scattering medium of 15 μm × 12 μm, more than half of the incident light can be delivered and absorbed at the target position. The enhancement in absorption is more than four times higher than that with random input. This value can be even higher for smaller absorption regions. We also quantify the effectiveness of the method and show that it is inversely proportional to the openness of the scattering medium. This result is potentially useful for targeted light delivery inside scattering media with local absorption.
© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement
Corrections
17 June 2021: Typographical corrections were made to Figs. 8–11.
1. Introduction
Optical scattering inherent in biological tissue is a long-standing problem that limits the operational range of optical imaging, manipulation, and therapy [1]. Recent development in wavefront shaping allows coherent control of scattered light using a spatial light modulator (SLM), making deep-tissue non-invasive optical imaging and optogenetics attainable [2–6]. Successful operation of wavefront shaping relies on the determination of transmission matrices of scattering media, through direct measurement [7–13], feedback-based approach [14–20], or optical phase conjugation [21–27]. Here, the transmission matrix is a mathematical model to connect the input and the output field. With the knowledge of the transmission matrix, one can synthesize arbitrary field patterns through the scattering media. In the special case where light is focused to a single speckle through scattering media, intensity enhancement over tens of thousands can be achieved [28–30]. For lossless scattering media, the eigenvalues and eigenchannels of transmission matrices can be described by the random matrix theory (RMT) [31, 32]. Due to energy conservation, the probability density function of the eigenvalues τ follows a bimodal distribution that peaks at τ = 1 and τ = 0, corresponding to open channels and closed channels. In the past decades, researchers spent numerous efforts to couple shaped wavefront to these open channels, to enhance transmittance [33–36]. The best performance achieved to date is about four times enhancement for total transmittance through a 27-μm-thick scattering layer filled with ZnO nanoparticles [37].
In practice, scattering media are generally not lossless. Absorption does exist in scattering media, altering the behaviors of scattered light [38, 39]. In particular, it has been shown that absorption can significantly alter the propagation of power flow, especially for the open channels [40, 41]. Moreover, uniform absorption of scattering media can turn the diffusive transport of light in the maximum transmission channel into quasi-ballistic [39]. Consequently, the distribution of transmission eigenvalues follows the quarter-circle law rather than bimodal [7, 34]. In many practical situations, targeted light deposition inside scattering media is essentially important, with the potential to benefit signal strength and imaging resolution for deep-tissue photoacoustic and fluorescence microscopy [3–5]. This desire is closely related to an effect named coherent enhanced absorption (CEA) [42]. Chong and Stone showed that for a terminated waveguide geometry with reflection only, CEA is closely related to the smallest reflection eigenvalue [42]. Later, it was shown that CEA could be made broadband due to long-range correlations, causing the spectral degrees of freedom to scale as the square root of the bandwidth [43]. Despite these accomplishments, the performance of CEA remains largely unexplored for scattering media with open boundaries. In this situation, one would expect to collect and time-reversed the scattered light from all possible directions as much as possible. Our study explores this situation by showing how to establish a numerical model assisted by the genetic algorithm (GA) [16] to achieve effective CEA, when the shaped wavefront is synthesized on only one side of the scattering medium with finite controls.
Although being the golden standard to numerically solve Maxwell’s equations of electromagnetic waves, the method of finite-difference time-domain (FDTD) was rarely employed as the numerical tool to study scattering media. This condition is due to the large computational resource required to simulate macroscopic scattering media with microscopic structures. To simplify the computational burden, many previous works employed either a tight-binding model [41] or a recursive Green’s function method [35, 40, 42, 43] numerically. To the best of the knowledge, only Choi’s group reported modeling the transmission matrix in 2011 [44] and the time-gated reflection matrix in 2018 [45] using the FDTD simulations. Besides, the method of finite-difference frequency-domain was employed to construct the transmission and reflection matrices for a disordered waveguide with perfectly reflecting sidewalls [39]. In this study, we developed a numerical framework for scattering media with open boundaries and construct both the transmission and reflection matrices through the FDTD. Singular value decomposition (SVD) is employed to analyze the properties of the transmission and reflection matrices with and without absorption. Then, we develop a GA-assisted numerical model to find the appropriate wavefront to realize CEA inside scattering media. Numerical results show that by treating the total collected transmittance and reflectance as the fitness function to be minimized, CEA can be realized for a single absorber buried inside the scattering medium with more than four times enhancement in absorption, when the shaped wavefront is synthesized at only one side. Furthermore, for three discrete absorbers buried inside the scattering media, we found that different realizations of CEA can be observed, exhibiting rich behaviors. We envision that this work can be beneficial to experimental approaches to efficiently deliver light to target absorptive positions.
2. Numerical model of transmission and reflection matrices
We start by describing how to construct transmission and reflection matrices. To solve Maxwell equations, Lumerical FDTD Solutions is employed. To ease computational complexity, in this study, we restrict ourselves in solving numerical solutions for two-dimensional geometry that lies in the x-y plane and is infinitely long along the z-axis. A schematic of the scattering medium we studied in this work is shown in Fig. 1(a), which consists of many dielectric cylinders randomly distributed in a vacuum as scatterers. By varying the refractive index n and the diameter d, we can control the scattering anisotropy g of the scattering. The computational region has a width W and a length L. The filling ratio, which represents the volume density, denotes the scattering mean free path lt.
Fig. 1. Numerical model of the computational environment. (a) Schematic of the scattering medium. Dielectric cylinders are randomly distributed in the computational region. (b) and (c) Amplitude and phase of the transmission matrix. (d) A typical illustration of the intensity enhancement, when the incident field is synthesized as the conjugation of one row of the transmission matrix. (e) and (f) Amplitude and phase of the reflection matrix.
In the first example, we choose n = 2.0 and d = 200 nm. The width and length of the scattering region are chosen as W = 15 μm and L = 12 μm. The filling factor is estimated as 50.3%. When the incident wavelength λ is 600 nm and the electric field is polarized along the z-axis (TM polarization), it can be theoretically estimated from Mie theory that g ≈ 0.6162 [46]. Following the produces described in Ref. [44], the scattering mean free path is numerically estimated to be lt ≈ 328 nm. Thus, the length of the scattering medium is about 37 times thicker than lt, which is sufficiently scattering to scramble the input wavefront. During simulations, perfectly matched layers (PML) are set for all boundaries to mimic open boundaries. The size of the mesh grid is set to 10 nm (λ/60), which guarantees computational accuracy.
To reconstruct transmission and reflection matrices, a series of the input field ${E_{\textrm{in}}}({{k_x} > 0,{k_y}} )$ are sent to probe the scattering medium. Here, ${k_0} = 2\pi /\lambda $ is the wave vector in a vacuum, while ${k_x} = {k_0}\cos \theta $ and ${k_y} = {k_0}\sin \theta $ are the corresponding components along the x- and y- directions, respectively. $\theta $ is the illumination angle of the input field. Throughout this work, we set the angle θ with a range from -73.7° to 73.7°, as larger angles could make difficulties in the convergence. Although this choice makes the constructed matrices incomplete and increases the openness of the scattering medium as well, the computational burden is greatly alleviated. Thus, the degree of freedom on the illumination side is $W/({\lambda /2 \times \textrm{sin}73.7^\circ } )\approx 48$. To make the illumination angle symmetric concerning θ = 0, 49 discrete incident angle θ are chosen with $\sin \theta $ being equally spaced. Then, we transverse the input light and record the corresponding scattered light at all boundaries. A built-in module of the software is applied to compute far-field projections ${E_{\textrm{out}}}({{k_x}^{\prime},{k_y}^{\prime}} )$ based on the uniqueness theorem of the electromagnetic waves. The transmission matrix in the space of the wave vector is constructed using the following relation
Similarly, the reflection matrix in the space of the wave vector is constructed using the following equation
3. Properties of transmission and reflection matrices
3.1 Lossless scattering media
Having described the numerical scheme to construct transmission and reflection matrices, we first examine the properties of these matrices. According to the RMT, the transmission matrix can be decomposed through SVD
where τ is a diagonal matrix. Each diagonal element of τ is a non-negative real number and the square of the i-th component $\tau {({i,i} )^2}$ physically represents the transmittance of the i-th transmission eigenchannel. Ut and Vt are unitary matrices. ${U_t}({k_y^{\prime},i} )$ and ${V_t}({i,{k_y}} )$ map the i-th transmission eigenchannel to the transmitted field with wave vector $k_y^{\prime}$ and the input field with wave vector ${k_y}$ to the i-th transmission eigenchannel, respectively. Based on the RMT, the density of the transmittance T for lossless scattering media follows the bimodal distribution [31, 32]: where T is the ensemble-averaged transmittance. This distribution has two peaks at T = 1 and T = 0, corresponding to “open channels” and “closed channels”, respectively. These two types of channels are the natural results of energy conservation as light can be either transmitted or reflected. By numerically decomposing the transmission matrix constructed above, Fig. 2(a) plots the transmittance T for different transmission eigenchannels using red dots. As a comparison, transmittance computed based on RMT with $T = 0.182$ (the black dashed line) is also plotted using a blue curve. As shown in Fig. 2(a), transmittance obtained through SVD and RMT matches quite well, confirming the validity of the numerical model. We note here that since the presented distribution is different from the one predicted using the theory of filtering random matrix [49], the scattering geometry we considered here is quasi-closed with matrix elements being largely correlated. The only noticeable discrepancy occurs around the 20th channel, possibly due to the energy leaked from the side of the scattering medium. Moreover, only several transmission eigenchannels exhibit considerable transmittance, while the transmittance for the rest is quite small. This observation indicates that a thorough investigation of the scattering medium requires knowledge of the reflection matrix as well, which is usually omitted in previous studies.Fig. 2. Properties of the eigenvalues and eigenchannels for transmission and reflection matrices for a lossless scattering medium. (a) Eigenvalue distribution of the transmission matrix is displayed in descending order. Results obtained through the method of the FDTD (red dots) and the RMT (blue curve) are compared. (b) Eigenvalue distribution of reflection matrix displayed in ascending order. (c) Correlation coefficients between each pair of the transmission and reflection eigenchannels with the same indices. (d) Correlation coefficients between the first transmission eigenchannel and every reflection eigenchannel.
Borrowing the procedures to decompose the transmission matrix, we can perform SVD to the reflection matrix as well
Similarly, $\rho $ is a diagonal matrix while Ur and Vr are unitary matrices. Each diagonal element of $\rho $ is a non-negative real number and the square of the i-th component $\rho {({i,i} )^2}$ physically represents the reflectance of the i-th reflection eigenchannel. ${U_r}({k_y^{\prime},i} )$ and ${V_r}({i,{k_y}} )$ map the i-th reflection eigenchannel to the reflected field with wave vector $k_y^{\prime}$ and the input field with wave vector ${k_y}$ to the i-th reflection eigenchannel, respectively. Figure 2(b) plots the reflectance R for different reflection eigenchannels using red dots. The displaying order of the reflection eigenchannels is reversed, showing a complementary structure to the one in Fig. 2(a). The average reflectance is 0.239 (the black dashed line). Since the summation of the averaged transmittance and reflection is far below 1, indicating that a considerable amount of energy is leaked from the side of the scattering medium. This condition is commonly seen for scattering media during experiments, where the collected transmitted light and reflected light occupy only a small portion of the total scattered light [50].Then, we explore the correlation coefficient between the transmission and reflection eigenchannels. Figure 2(c) illustrates the correlation coefficient, defined as the absolute value of the 1st order field correlation, between the transmission and reflection eigenchannel with the same indices. The highest correlation coefficient 0.792 is achieved between the 1st transmission and reflection eigenchannels, corresponding to the T ≈ 1 and R ≈ 0, respectively. In contrast, the correlation coefficient between the last transmission and reflection eigenchannels, corresponding to the T ≈ 0 and R ≈ 1, decreases to 0.702. This asymmetric result between transmission and reflection is possibly due to the surface reflection that is only dominant on the reflection side. Nonetheless, high correlation coefficients can always be achieved for all pairs, with an average value of 0.503. This result is expected, as a strong(weak) transmission generally leads to a weak(strong) reflection. We also compute the correlation coefficients between the 1st transmission eigenchannel and every reflection eigenchannel. As shown in the figure, all reflection eigenchannels other than the first show relatively small correlation coefficients, with an average value of only 0.153. These observations are consistent with the experimental results reported in Ref. [36].
To visualize the power flow inside the scattering medium for different eigenchannels, we generate two-dimensional vector plots of Poynting vectors. Figures 3(a) and (b) show the power flow when synthesizing the incident field to excite the 1st transmission and the last reflection eigenchannels, exhibiting significant differences in terms of transmittance and reflectance. Specifically, for Fig. 3(a), the transmittance, reflectance, and the leaked light occupy 93.3%, 2.1%, and 4.6% of the input energy, respectively. In contrast, these numbers are 2.2%, 96.1%, 1.7% for Fig. 3(b). We note that as the openness of the scattering medium increases, behaviors to achieve near-unity transmission or reflection will become more difficult due to increased leaked light. As a control, Fig. 3(c) shows the power flow for a random incident wavefront, exhibiting comparative transmission and reflection on both sides of the computational region.
Fig. 3. The power flow inside the scattering media, denoted by two-dimensional vector plots of Poynting vectors. (a) The incident wavefront is synthesized to excite the 1st transmission eigenchannel (T ≈ 1, “open channel”). (b) The incident wavefront is synthesized to excite the last reflection eigenchannel (R ≈ 1, “closed channel”). (c) A random wavefront is synthesized and sent into the scattering medium.
3.2 Scattering media with absorption
We then investigate the situation when a local absorber is introduced inside the scattering medium, which is commonly seen in biological tissue. Numerically, the local absorber can be simulated by adding a nonzero imaginary part to the refractive index. As shown in Fig. 4(a), we first consider the case of a single absorber with a diameter of D = 5 μm. As a proof of concept, we artificially set the refractive index of the absorbing region as n = 2 + i8. Figures 4(b) – (e) show the amplitude and phase of the constructed transmission and reflection matrices. Compared to the lossless case, one could not visually tell the statistical difference after introducing this local absorber.
Fig. 4. (a) Schematics of the scattering medium with a single absorber. (b) and (c) Amplitude and phase of the transmission matrix. (d) and (e) Amplitude and phase of the reflection matrix.
Then, we employ SVD to decompose these matrices to investigate how this introduced absorber alters the properties of the eigenchannels. Figure 5(a) plots the computed transmittance of eigenchannels using red dots. Due to absorption, the averaged transmittance decreases to 0.107 (the black dashed line). As shown in the figure, the existence of the local absorber significantly reduced the transmittance for all the eigenchannels. Even for the 1st eigenchannel, the transmittance drops to 0.842. Moreover, the transmittance is considerably smaller than that of its lossless companion computed using the RMT. Similarly, Fig. 5(b) plots the computed reflectance of the eigenchannels using red dots, and the averaged reflectance also decreases to 0.131 (the black dashed line) due to absorption. Figure 5(c) illustrates the correlation coefficient between the transmission and reflection eigenchannels with the same indices. Different from the lossless case, only the 1st transmission and reflection eigenchannels are strongly correlated, while the correlation coefficients for the rest pairs significantly drop. As a result, the mean correlation coefficient drops to 0.206. Figure 5(d) plots correlation coefficients between the 1st transmission eigenchannel and every reflection eigenchannel, exhibiting similar behaviors as that of Fig. 2(d).
Fig. 5. Properties of the eigenvalues and eigenchannels for transmission and reflection matrices for a single absorber buried inside a scattering medium. (a) Eigenvalue distribution of the transmission matrix is displayed in descending order. Results obtained through the method of the FDTD (red dots) and the RMT (blue curve) are compared. (b) Eigenvalue distribution of reflection matrix displayed in ascending order. (c) Correlation coefficients between each pair of the transmission and reflection eigenchannels with the same indices. (d) Correlation coefficients between the first transmission eigenchannel and every reflection eigenchannel.
4. Coherent enhanced absorption in scattering media
4.1 Single absorber
In many practical situations, there is a strong need to deliver energy to a target position inside scattering media. Based on the time-reversal principle, coherent perfect absorption, which is an ideal situation, requires the incident light to be modulated with phase and amplitude from all directions. However, this situation is not possible in practice. In this section, we seek the optimum solution for the incident field that is generated on only one side of the scattering medium, with the purpose to enhance local absorption. With the knowledge of the transmission and reflection matrices, the GA [16] was employed in assisting in searching for the desired wavefront.
The flowchart of the GA is schematically shown in Fig. 6. The GA begins by generating an initial population of NP phase masks. The fitness of each phase mask is evaluated as the sum of the total transmittance and reflectance and the goal is to minimize the fitness. The underlying assumption is that light is likely to be absorbed when both the transmittance and reflectance are small. To breed offspring, two parents, pa and ma, are randomly selected from the population. A new offspring is produced as pa · T + ma · (1 - T), where T is a random binary template generated with a parental mask ratio pc. In the following, an adaptive mutation operation is performed to avoid being trapped into the local optimum. After mutation, the fitness of the new offspring is measured. The steps of generating new offspring are repeated many times until a generational condition is fulfilled. A fixed number of 100 generations are performed for the GA in this work. The population of each generation is fixed at NP = 80. For the breeding operation, the mask ratio pc is 0.5. For the mutation operation, the initial mutation rate, the final mutation rate, and the decay factor are set to be 0.1, 0.0025, and 200, respectively.
Figure 7(a) shows a typical running process of the GA and plots the T + R as a function of the number of measurements using a blue curve. In the beginning, the fitness is 1 and its corresponding power flow is shown in Fig. 7(b), exhibiting randomly distributed power flow. By continuously running the GA, T + R keeps decreasing and drops to around 0.5 after 2000 measurements. The corresponding power flow is shown in Fig. 7(c), showing an enhanced energy density within the absorptive region. The T + R is further reduced with an increased number of generations. At roughly 4900 measurements, the corresponding power flow is shown in Fig. 7(d), showing a further increased energy density within the absorption region. The T + R measured through FDTD simulation when synthesizing the corresponding wavefront is denoted as a red dot, which is close to the one computed using the GA and the matrices. Before proceeding, we note that a small T + R does not necessarily mean an efficient CEA at the target position, as light can be leaked from the sides as well, as shown in Fig. 6(e). This situation is commonly seen as a local optimum in the GA, indicating the necessity to check the distribution of power flow. In practice, this situation can be identified and avoided by placing several optical detectors at the sides of the scattering medium.
Fig. 7. Results of the GA-assisted CEA for the scattering medium with a single absorber. (a) A typical running process of the GA, showing T + R as a function of the number of measurements. (b) Power flow of a random input wavefront. (c) Power flow of the weak absorbing case. (d) Power flow of the strong absorbing case. (e) Power flow of the case with light leaking from the sides of the scattering medium.
For better comparison, the transmittance, reflectance, the percentage of the leaked energy from sides, and the absorption for several representative cases corresponding to Figs. 7(b)-(e) are quantified and listed in Tab. 1. The absorption in case (c) is more than twice that in case (b), while the absorption in case (d) is roughly four times that in case (b). Remarkably, more than half of the incident light is absorbed in case (d), indicating CEA inside the scattering medium can be effectively achieved by manipulating the wavefront of incident light from only one side.
Table 1. Quantitative analysis of CEA for a single absorber
Without setting internal monitors, it seems challenging to differentiate the cases of CEA and leakage. To address this issue, we first consider in what condition, the incident light is likely to be absorbed by the cylinder or leaked from the sides. Figure 8 (a) plots the percentages of light being leaked or absorbed for light with different incident angles. Consistent with physical intuition, light with small incident angles is inclined to be absorbed, while light with large angles is likely to leak from the sides. Here, we set arctan(W/2L) ≈ 36.9° as the demarcation point, by considering the incident light at the center. Therefore, we speculate to use this observation as a criterion to differentiate the cases of CEA and leakage, by examining the composition of the synthesized optical field. Specifically, we define a(θ) as the absolute value of the coefficient that corresponds to incident light with an angle θ. For the scattering medium considered here, we stipulate that if the largest three a(θ) are within the range of -36.9° to 36.9°, the optimized field probably leads to CEA. In contrast, if all the largest three a(θ) are beyond this range, the optimized field probably corresponds to leakage. The reason we chose three values is that for CEA, we expect the largest one locates around θ ≈ 0° while the second largest two are located at both sides near the largest one. To validate this assumption, Figs. 8(b)-(d) plot the corresponding amplitudes of the coefficients a(θ) for the cases of weak, strong, and leakage in Fig. 7, respectively. For the case of weak absorption, the largest three values correspond to incident angles of -26.1°, -13.9°, and 49.5° (Fig. 8(b)). Since not all three angles are within or out of the range of -36.9° to 36.9°, this synthesized field corresponds to neither CEA nor leakage. This conclusion is consistent with the fact that it is indeed a local optimum obtained through the GA. For the case of strong absorption, the largest three values correspond to incident angles of -9.2°, 0°, and 21.1°. As a comparison, for the case of leakage, the largest three values correspond to incident angles of -66.9°, 39.8°, and 61.6°. These results confirm that examining incident angles of the largest three amplitudes of coefficients is a good criterion to differentiate different cases without employing internal monitors.
Fig. 8. (a) Percentages of light being leaked or absorbed for input light with different incident angles. (b)-(d) Plots of the corresponding amplitudes of the coefficients a(θ) for cases of weak, strong, and leakage in Fig. 7, respectively.
Besides, a systematic study on the statistical distribution of a(θ) can also be used to differentiate the cases of CEA and leakage. First, we modified a(θ) by eliminating the components whose values are below the mean and only keeping those with large values. Then, we define a step function s(θ), whose value is set to 1 if $|\theta |\le 36.9^\circ $ and is set to 0 otherwise. The classification can be made by computing the correlation between modified a(θ) and s(θ). If modified a(θ) and s(θ) are positively correlated, we can conclude that the synthesized field is likely to be classified as CEA, as a(θ) concentrates more on small angles. Otherwise, for negative correlation, the synthesized field likely leads to the case of leakage, as a(θ) concentrates more on large angles. Numerically, we found the correlations for cases shown in Figs. 8(c) and (d) are 0.438 and -0.707, respectively, confirming the above statement.
4.2 Three distinct absorbers
We further consider the situation where three distinct absorbers are located within the scattering medium. Figure 9 shows the schematic of the computational region, where three absorptive regions have the same diameter of 5 μm and the same refractive index of 2 + i8. Following the same procedure, the transmission and reflection matrices are constructed with amplitude and phase shown in Figs. 9(b)-(g).
Fig. 9. (a) Schematics of the scattering medium with three distinct absorbers. (b) and (c) Amplitude and phase of the transmission matrix. (d) and (e) Amplitude and phase of the reflection matrix.
By employing SVD, the transmittance and reflectance of the corresponding eigenchannels are plotted in Figs. 10 (a) and (b). Compared to the case with only one absorber, the averaged transmittance and reflectance become smaller. Correlation coefficients between different eigenchannels are also shown in Figs. 9(c) and (d), exhibiting similar statistical behaviors like that of the single absorber case.
Fig. 10. Properties of the eigenvalues and eigenchannels for transmission and reflection matrices for three distinct absorbers buried inside a scattering medium. (a) Eigenvalue distribution of the transmission matrix is displayed in descending order. Results obtained through the method of the FDTD (red dots) and the RMT (blue curve) are compared. (b) Eigenvalue distribution of reflection matrix displayed in ascending order. (c) Correlation coefficients between each pair of the transmission and reflection eigenchannels with the same indices. (d) Correlation coefficients between the first transmission eigenchannel and every reflection eigenchannel.
We also apply the GA to achieve CEA in this geometry by manipulating the incident wavefront from only one side. The parameters of the GA are set to be the same as the single-absorber case. Interestingly, several different conditions can be realized, with corresponding power flows shown in Figs. 11(a)-(d). Figure 11(a) corresponds to the case with a randomly generated incident wavefront, and no significant difference in the energy density in terms of position is observed. Figure 11(b) shows the case that light is mostly absorbed by only one of the absorbers, while Fig. 11(c) shows the case that light is mostly absorbed by two absorbers. Moreover, light can also be leaked out from the sides of the scattering medium, as shown in Fig. 11(d). Quantitative analysis of the energy distribution for each case is detailed in Tab. 2. We found that due to the increased area of the absorptive region, the enhancement in the absorption is worse compared to that in the single absorber case when the number of independent controls remains the same. Moreover, although both cases (b) and (c) realized about half of the incident energy being absorbed, the absorption in (b) is larger than that in (c). These results indicate that CEA is more challenging with the enlarged absorption area. Nonetheless, several realizations of CEA inside scattering media can be effectively achieved by manipulating incident wavefront from only one side.
Fig. 11. Results of the GA-assisted CEA for a scattering medium with three distinct absorbers. (a) Power flow of a random input wavefront. (b) Power flow of the case with light being absorbed by only one of the absorbers. (c) Power flow of the case with light being absorbed by two absorbers. (d) Power flow of the case with light leaking from the sides of the scattering medium.
Table 2. Quantitative analysis of CEA for three distinct absorbers
5. Discussion
5.1 Diameters of the absorber
For the single absorber case, we further discuss how the areas of the absorption region affect the performance of CEA. To examine this issue, we perform the above procedures to realize CEA for a single absorber with diameters of 2 and 8 μm. The results are summarized and plotted in Fig. 12. As shown using red dots, the absorption realized through GA-assisted CA also increases with an increased diameter. However, since the absorption for unoptimized input light also increases, a fair way to describe the effectiveness of GA-assisted CEA is to introduce a new parameter η. This parameter is defined as the energy ratio between the absorption of GA-assisted CEA and the ensemble-averaged absorption with random input. In contrast to the absolute absorption, η decreases considerably with the increased diameter (blue squares). This observation again confirms that it is generally more difficult to manipulate a larger target for the same number of input controls using wavefront shaping. Nonetheless, we note that η does not scale down with the inverse square of the diameter, showing the promise to realize CEA for a large region with finite controls. Such a capability is probably due to the long-range correlations in the scattering region [43].
Fig. 12. The performance of CEA for a single absorber with different diameters. The red dots and blue squares represent the absorption with GA-assisted CEA and the enhancement η, respectively.
5.2 Openness of scattering media
Another important factor to be investigated is the openness of the scattering medium, as the extent of the leaked light can significantly deteriorate the performance of CEA. To quantify this effect, we considered several scattering media with different thicknesses. Figure 13(a) plots the simulation results for scattering media with thicknesses of 12, 20, 25, and 30 μm. the ensemble-averaged leakage for unoptimized input light, which can be used to represent openness, are about 27.0%, 34.1%, 38.4%, and 41.8%, respectively, which are plotted using red dots. Due to the huge computational burden of the FDTD simulation, which will be shown in the next subsection, we did not further study thicker scattering media. For these scattering media, η of the GA-assisted CEA are plotted using blue squares in Fig. 13(a), which decays as the thickness increases. Figure 13(b) directly shows how η evolves as a function of leakage. Notably, the effectiveness of this method decreases as the openness of the scattering medium increases, indicating our method only works for systems with finite leakage. Thus, to realize CEA for practical experimental conditions with open geometries, effective approaches such as employing internal guide stars need to be explored in the future.
Fig. 13. The performance of CEA for a single absorber buried within scattering media with different thicknesses. (a) The red dots and blue squares represent the ensemble-averaged leakage with an unoptimized input field and the enhancement η of the GA-assisted CEA, respectively. (b) Enhancement η as a function of leakage.
5.3 Computational time cost for the FDTD and the GA
Even using commercial software to compute such simplified two-dimensional geometry, the FDTD computation is extremely time-consuming to guarantee sufficient accuracy. Here, we summary the typical computational time cost to realize the GA-assisted CEA in Tab. 3. The computer we used for simulation has an i7-10700 @ 2.90GHz CPU and a Kingston DDR4 2666 MHz 64GB RAM. As shown in Tab. 3, even for the scattering medium with a thickness of 12 μm, it took roughly 90 hours. Moreover, the computational time increases nonlinearly with the increased thickness and reaches 320 hours to simulate the scattering medium with a thickness of 30 μm. Compared to the FDTD simulation, the computational time cost for the GA to search for the optimum solution is much less, which is roughly 20 seconds regardless of the size of the scattering medium.
Table 3. Computational time cost for different scattering media
6 Conclusion
In conclusion, we develop a numerical scheme to construct transmission and reflection matrices of scattering media with open boundaries. Based on these matrices, the properties of the eigenchannels and eigenvalues of scattering media with and without absorption are investigated. To realize effective CEA for scattering media with different types of absorption, the GA is employed to minimize the total transmittance and reflectance simultaneously, allowing the synthesis of the incident wavefront from only one side of the scattering medium. We further introduce inhomogeneous absorption to our simulation. In the case of one absorber, intensity enhancement up to four times can be achieved, resulting in the absorption of half of the incident energy. In the case of three absorbers, different realizations such as enhancement in one of the absorbers and two of the absorbers have been demonstrated. The effectiveness of the GA-assisted CEA is further investigated for absorbers with different diameters and scattering media with different openness. This work is anticipated to pave a way for future experimental approaches to achieve targeted energy deposition.
Funding
National Key Research and Development Program of China (2018YFB1802300); National Natural Science Foundation of China (12004446); Fundamental and Applied Basic Research Project of Guangzhou (202102020603).
Disclosures
The authors declare no conflicts of interest.
Data availability
Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.
References
1. V. Ntziachristos, J. Ripoll, L. V. Wang, and R. Weissleder, “Looking and listening to light: the evolution of whole-body photonic imaging,” Nat. Biotechnol. 23(3), 313–320 (2005). [CrossRef]
2. S. Popoff, G. Lerosey, M. Fink, A. C. Boccara, and S. Gigan, “Image transmission through an opaque material,” Nat. Commun. 1(1), 81 (2010). [CrossRef]
3. 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(1), 928 (2012). [CrossRef]
4. K. Si, R. Fiolka, and M. Cui, “Fluorescence imaging beyond the ballistic regime by ultrasound-pulse-guided digital phase conjugation,” Nat. Photonics 6(10), 657–661 (2012). [CrossRef]
5. D. B. Conkey, A. M. Caravaca-Aguirre, J. D. Dove, H. Ju, T. W. Murray, and R. Piestun, “Super-resolution photoacoustic imaging through a scattering wall,” Nat. Commun. 6(1), 7902 (2015). [CrossRef]
6. H. Ruan, J. Brake, J. E. Robinson, Y. Liu, M. Jang, C. Xiao, C. Zhou, V. Gradinaru, and C. Yang, “Deep tissue optical focusing and optogenetic modulation with time-reversed ultrasonically encoded light,” Sci. Adv. 3(12), eaao5520 (2017). [CrossRef]
7. S. M. Popoff, G. Lerosey, R. Carminati, M. Fink, A. C. Boccara, and S. Gigan, “Measuring the Transmission Matrix in Optics: An Approach to the Study and Control of Light Propagation in Disordered Media,” Phys. Rev. Lett. 104(10), 100601 (2010). [CrossRef]
8. T. Chaigne, O. Katz, A. C. Boccara, M. Fink, E. Bossy, and S. Gigan, “Controlling light in scattering media non-invasively using the photoacoustic transmission matrix,” Nat. Photonics 8(1), 58–64 (2014). [CrossRef]
9. M. Kim, W. Choi, Y. Choi, C. Yoon, and W. Choi, “Transmission matrix of a scattering medium and its applications in biophotonics,” Opt. Express 23(10), 12648–12668 (2015). [CrossRef]
10. J. Yoon, K. Lee, J. Park, and Y. Park, “Measuring optical transmission matrices by wavefront shaping,” Opt. Express 23(8), 10158–10167 (2015). [CrossRef]
11. J. Xu, H. Ruan, Y. Liu, H. Zhou, and C. Yang, “Focusing light through scattering media by transmission matrix inversion,” Opt. Express 25(22), 27234–27246 (2017). [CrossRef]
12. T. Zhao, L. Deng, W. Wang, D. S. Elson, and L. Su, “Bayes’ theorem-based binary algorithm for fast reference-less calibration of a multimode fiber,” Opt. Express 26(16), 20368–20378 (2018). [CrossRef]
13. X. Wei, Y. Shen, J. C. Jing, A. S. Hemphill, Z. Xu, Z. Yang, Z. Yang, and L. V. Wang, “Real-time frequency-encoded spatiotemporal focusing through scattering media using a programmable 2D ultrafine optical frequency comb,” Sci. Adv. 6(8), eaay1192 (2020). [CrossRef]
14. I. M. Vellekoop and A. P. Mosk, “Focusing coherent light through opaque strongly scattering media,” Opt. Lett. 32(16), 2309–2311 (2007). [CrossRef]
15. O. Katz, E. Small, Y. Bromberg, and Y. Silberberg, “Focusing and compression of ultrashort pulses through scattering media,” Nat. Photonics 5(6), 372–377 (2011). [CrossRef]
16. D. B. Conkey, A. N. Brown, A. M. Caravaca-Aguirre, and R. Piestun, “Genetic algorithm optimization for focusing through turbid media in noisy environments,” Opt. Express 20(5), 4840–4849 (2012). [CrossRef]
17. P. Lai, L. Wang, J. W. Tay, and L. V. Wang, “Photoacoustically guided wavefront shaping for enhanced optical focusing in scattering media,” Nat. Photonics 9(2), 126–132 (2015). [CrossRef]
18. Z. Wu, J. Luo, Y. Feng, X. Guo, Y. Shen, and Z. Li, “Controlling 1550-nm light through a multimode fiber using a Hadamard encoding algorithm,” Opt. Express 27(4), 5570–5580 (2019). [CrossRef]
19. Q. Feng, F. Yang, X. Xu, B. Zhang, Y. Ding, and Q. Liu, “Multi-objective optimization genetic algorithm for multi-point light focusing in wavefront shaping,” Opt. Express 27(25), 36459–36473 (2019). [CrossRef]
20. D. Wu, L. Qin, J. Luo, Z. Liao, X. Chen, H.-C. Chui, and Y. Shen, “Delivering targeted color light through a multimode fiber by field synthesis,” Opt. Express 28(13), 19700–19710 (2020). [CrossRef]
21. M. Cui and C. Yang, “Implementation of a digital optical phase conjugation system and its application to study the robustness of turbidity suppression by phase conjugation,” Opt. Express 18(4), 3444–3455 (2010). [CrossRef]
22. T. R. Hillman, T. Yamauchi, W. Choi, R. R. Dasari, M. S. Feld, Y. Park, and Z. Yaqoob, “Digital optical phase conjugation for delivering two-dimensional images through turbid media,” Sci. Rep. 3(1), 1909 (2013). [CrossRef]
23. B. Judkewitz, Y. M. Wang, R. Horstmeyer, A. Mathy, and C. Yang, “Speckle-scale focusing in the diffusive regime with time reversal of variance-encoded light (TROVE),” Nat. Photonics 7(4), 300–305 (2013). [CrossRef]
24. D. Wang, E. H. Zhou, J. Brake, H. Ruan, M. Jang, and C. Yang, “Focusing through dynamic tissue with millisecond digital optical phase conjugation,” Optica 2(8), 728–735 (2015). [CrossRef]
25. H. Ruan, M. Jang, and C. Yang, “Optical focusing inside scattering media with time-reversed ultrasound microbubble encoded light,” Nat. Commun. 6(1), 8968 (2015). [CrossRef]
26. Y. Liu, C. Ma, Y. Shen, J. Shi, and L. V. Wang, “Focusing light inside dynamic scattering media with millisecond digital optical phase conjugation,” Optica 4(2), 280–288 (2017). [CrossRef]
27. C. Ma, J. Di, Y. Zhang, P. Li, F. Xiao, K. Liu, X. Bai, and J. Zhao, “Reconstruction of structured laser beams through a multimode fiber based on digital optical phase conjugation,” Opt. Lett. 43(14), 3333–3336 (2018). [CrossRef]
28. M. Jang, H. Ruan, H. Zhou, B. Judkewitz, and C. Yang, “Method for auto-alignment of digital optical phase conjugation systems based on digital propagation,” Opt. Express 22(12), 14054–14071 (2014). [CrossRef]
29. Y. Shen, Y. Liu, C. Ma, and L. V. Wang, “Sub-Nyquist sampling boosts targeted light transport through opaque scattering media,” Optica 4(1), 97–102 (2017). [CrossRef]
30. H. Yu, K. Lee, and Y. Park, “Ultrahigh enhancement of light focusing through disordered media controlled by mega-pixel modes,” Opt. Express 25(7), 8036–8047 (2017). [CrossRef]
31. O. N. Dorokhov, “On the coexistence of localized and extended electronic states in the metallic phase,” Solid State Commun. 51(6), 381–384 (1984). [CrossRef]
32. C. W. J. Beenakker, “Random-matrix theory of quantum transport,” Rev. Mod. Phys. 69(3), 731–808 (1997). [CrossRef]
33. I. M. Vellekoop and A. P. Mosk, “Universal Optimal Transmission of Light Through Disordered Materials,” Phys. Rev. Lett. 101(12), 120601 (2008). [CrossRef]
34. H. Yu, T. R. Hillman, W. Choi, J. O. Lee, M. S. Feld, R. R. Dasari, and Y. Park, “Measuring Large Optical Transmission Matrices of Disordered Media,” Phys. Rev. Lett. 111(15), 153902 (2013). [CrossRef]
35. S. M. Popoff, A. Goetschy, S. F. Liew, A. D. Stone, and H. Cao, “Coherent Control of Total Transmission of Light through Disordered Media,” Phys. Rev. Lett. 112(13), 133903 (2014). [CrossRef]
36. M. Kim, W. Choi, C. Yoon, G. H. Kim, S.-h. Kim, G.-R. Yi, Q. H. Park, and W. Choi, “Exploring anti-reflection modes in disordered media,” Opt. Express 23(10), 12740–12749 (2015). [CrossRef]
37. M. Kim, Y. Choi, C. Yoon, W. Choi, J. Kim, Q.-H. Park, and W. Choi, “Maximal energy transport through disordered media with the implementation of transmission eigenchannels,” Nat. Photonics 6(9), 581–585 (2012). [CrossRef]
38. A. G. Yamilov and B. Payne, “Interplay between localization and absorption in disordered waveguides,” Opt. Express 21(10), 11688–11697 (2013). [CrossRef]
39. S. F. Liew, S. M. Popoff, A. P. Mosk, W. L. Vos, and H. Cao, “Transmission channels for light in absorbing random media: From diffusive to ballistic-like transport,” Phys. Rev. B 89(22), 224202 (2014). [CrossRef]
40. S. F. Liew and H. Cao, “Modification of light transmission channels by inhomogeneous absorption in random media,” Opt. Express 23(9), 11043–11053 (2015). [CrossRef]
41. A. Yamilov, S. Petrenko, R. Sarma, and H. Cao, “Shape dependence of transmission, reflection, and absorption eigenvalue densities in disordered waveguides with dissipation,” Phys. Rev. B 93(10), 100201 (2016). [CrossRef]
42. Y. Chong and A. D. Stone, “Hidden black: Coherent enhancement of absorption in strongly scattering media,” Phys. Rev. Lett. 107(16), 163901 (2011). [CrossRef]
43. C. W. Hsu, A. Goetschy, Y. Bromberg, A. D. Stone, and H. Cao, “Broadband Coherent Enhancement of Transmission and Absorption in Disordered Media,” Phys. Rev. Lett. 115(22), 223901 (2015). [CrossRef]
44. W. Choi, A. P. Mosk, Q. H. Park, and W. Choi, “Transmission eigenchannels in a disordered medium,” Phys. Rev. B 83(13), 134207 (2011). [CrossRef]
45. S. Jeong, Y.-R. Lee, W. Choi, S. Kang, J. H. Hong, J.-S. Park, Y.-S. Lim, H.-G. Park, and W. Choi, “Focusing of light energy inside a scattering medium by controlling the time-gated multiple light scattering,” Nat. Photonics 12(5), 277–283 (2018). [CrossRef]
46. C. F. Bohren and D. R. Huffman, Absorption and Scattering of Light by Small Particles (Wiley:1983).
47. G. Huang, D. Wu, J. Luo, Y. Huang, and Y. Shen, “Retrieving the optical transmission matrix of a multimode fiber using the extended Kalman filter,” Opt. Express 28(7), 9487–9500 (2020). [CrossRef]
48. G. Huang, D. Wu, J. Luo, L. Lu, F. Li, Y. Shen, and Z. Li, “Generalizing the Gerchberg-Saxton algorithm for retrieving complex optical transmission matrices,” Photon. Res. 9(1), 34–42 (2021). [CrossRef]
49. A. Goetschy and A. D. Stone, “Filtering Random Matrices: The Effect of Incomplete Channel Control in Multiple Scattering,” Phys. Rev. Lett. 111(6), 063901 (2013). [CrossRef]
50. Z. Yaqoob, D. Psaltis, M. S. Feld, and C. Yang, “Optical phase conjugation for turbidity suppression in biological samples,” Nature Photonics 2(2), 110–115 (2008). [CrossRef]
