Naturally occurring layered mineral franckeite with anisotropic Raman scattering and third-harmonic generation responses

Vertically stacked van der Waals (vdW) heterostructures have introduced a unique way to engineer optical and electronic responses in multifunctional photonic and quantum devices. However, the technical challenges associated with the artificially fabricated vertical heterostructures have emerged as a bottleneck to restrict their proficient utilization, which emphasizes the necessity of exploring naturally occurring vdW heterostructures. As one type of naturally occurring vdW heterostructures, franckeite has recently attracted significant interest in optoelectronic applications, but the understanding of light–matter interactions in such layered mineral is still very limited especially in the nonlinear optical regime. Herein, the anisotropic Raman scattering and third-harmonic generation (THG) from mechanically exfoliated franckeite thin flakes are investigated. The observed highly anisotropic Raman modes and THG emission patterns originate from the low-symmetry crystal structure of franckeite induced by the lattice incommensurability between two constituent stacked layers. The thickness-dependent anisotropic THG response is further analyzed to retrieve the third-order nonlinear susceptibility for franckeite crystal. The results discussed herein not only provide new insights in engineering the nonlinear light–matter interactions in natural vdW heterostructures, but also develop a testbed for designing future miniaturized quantum photonics devices and circuits based on such heterostructures.

shown as electrocatalysts for energy related applications to improve hydrogen evolution performance 29 . The optical properties of layered franckeite have also been studied including the refractive index in the visible spectrum determined by a quantitative analysis of the optical contrast spectra, as well as the third-order nonlinear optical properties characterized by Z-scan and spatial self-phase modulation techniques 23,30 . In addition, it has been recently shown that franckeite exhibits the spontaneous symmetry breakdown due to the PbS-like and SnS 2 -like lattice incommensurability, which results in anisotropic electrical, vibrational and optical responses 24 . Although many research efforts have been devoted to studying the nonlinear optical properties of 2D layered vdW materials 31 , there is no comprehensive study presented on the anisotropic nonlinear optical responses in the incommensurate heterostructures of franckeite yet.
Motivated by this, herein we demonstrate the anisotropic Raman scattering response as well as the anisotropic third-harmonic generation (THG) response from mechanically exfoliated franckeite thin flakes. The angleresolved polarized Raman spectroscopy and the polarization-dependent THG emission are utilized in probing the anisotropic natural vdW heterostructures and identifying the rippling direction of franckeite crystals. The Raman modes in franckeite crystal are found highly anisotropic under both parallel and perpendicular polarization configurations. It is observed that the THG emission pattern is also highly anisotropic with respect to the incident linearly polarized pump beam, arising from the symmetry-broken periodic ripples in the crystal structure of franckeite induced by the lattice incommensurability and in-plane spatial modulation of vdW interaction between the consecutively stacked PbS-like and SnS 2 -like layers. The measured anisotropic THG response is further corroborated with the theoretical nonlinear susceptibility model, which enables the determination of the relative magnitudes of the third-order nonlinear susceptibility tensor elements for franckeite crystal. The effect of the franckeite flake thickness on the THG emission power is further explored to determine the value of the third-order nonlinear susceptibility. We anticipate that these results will provide new insights into the comprehensive understanding of light-matter interactions in natural vdW heterostructures and in realizing advanced quantum photonics devices for future applications in integrated photonic circuits, polarization-based entangled photon generation, encrypted optical signal processing, and quantum information science.

Results
Anisotropic Raman response of franckeite thin flakes and identification of crystal rippling direction. Franckeite crystal consists of the repetitive stacking of two vdW layers along the c-axis, one is the PbS-like layer (Q-layer) and the other is the SnS 2 -like layer (H-layer), interacting together with vdW forces 21,22 . Figure 1a illustrates the simplified version of atomic layer arrangement (Q-layer and H-layer) in a typical franckeite crystal. The Q-layer consists of four atomic layers of sulfide compound with the general molecular formulation as MX, where M = Pb 2+ , Sn 2+ or Sb 3+ and X = S with the unit cell parameters a = 5.84 Å, b = 5.90 Å, and c = 17.3 Å 32 . The H-layer is composed of the octahedrons of sulfide compound with the general molecular formulation of MX 2 , where M = Sn 4+ or Fe 2+ and X = S with the unit cell parameters a = 3.68 Å, b = 6.32 Å, and c = 17.3 Å 32 . Franckeite has a triclinic crystal structure belonging to the P1 space group. The lattice incommensurability of two consecutively stacked PbS-like and SnS 2 -like layers introduces a strain-induced structural deformation, which develops symmetry-broken periodic ripples with a period of around 4.7 nm along the b-axis 24,32 . Such one-dimensional rippling of the heterostructure will give rise to the in-plane anisotropic physical properties for franckeite crystal in the a-b plane. Franckeite thin flakes are prepared by the typical mechanical exfoliation process using Nitto tape (SPV 224) from bulk natural franckeite mineral (from San José mine, Oruro, Bolivia) and transfer to a glass substrate. Figure 1b Fig. 1c,e, where the flake thickness is found approximately 21 nm and 138 nm, respectively.
The crystal axes of franckeite crystal can be identified from the Raman spectroscopy. First, the Raman spectra from the 138 nm-thick franckeite flake are characterized and analyzed. Figure 1f shows the measured Raman spectra with a 632.8 nm He-Ne laser as the excitation source in three different polarization analyzer configurations. The total Raman spectrum is measured without the analyzer after the sample. In the parallel or perpendicular configuration, the analyzer direction is set to be parallel or perpendicular to the linear polarization direction of the excitation beam. The Raman spectrum of franckeite exhibits several distinct Raman peaks within 70-400 cm − 3 38,39 . The 78 cm −1 peak is assigned as the A g mode of Sb 2 S 3 (73 cm −1 ). The peak at 93 cm −1 is attributed to the A g mode of SnS (95 cm −1 ). The 140 cm −1 peak corresponds to a combination of phonon modes of both SnS 2 lattices (the second-order effects at 140 cm −1 ) and PbS lattices (transverse acoustic and optical phonon modes at 154 cm −1 ). The 183 cm −1 peak is related to a combination of the A g mode of SnS (192 cm −1 ) and the A g mode of Sb 2 S 3 (191 cm −1 ). The 200 cm −1 peak is attributed to a combination of the E g mode of SnS 2 (202 cm −1 ) and the longitudinal optical phonon modes of PbS lattices (204 cm −1 ). The 258 cm −1 peak represents the B 1g /B 3g mode of Sb 2 S 3 (238 cm −1 ), while the 284 cm −1 peak belongs to the A g mode of Sb 2 S 3 (283 cm −1 ). The peak at 324 cm −1 is assigned as a combination of the A g mode of Sb 2 S 3 (312 cm −1 ) and the A 1g mode of SnS 2 (315 cm −1 ). Finally, two weak Raman peaks at 350 cm −1 and 368 cm −1 belong to the E g mode (337 cm −1 ) and the A g mode (372 cm −1 ) of FeS 2 , respectively. In addition, it is noticed that the majority of modes (93, 140, 183, 200, 259, 284 and 324 cm −1 ) appear in both the parallel and perpendicular configurations, and only few modes (78, 350 and 368 cm −1 ) can be observed in either configuration. The effect of the flake thickness variation on the Raman shift is also probed. Figure 1g  www.nature.com/scientificreports/ configuration. No significant change in the Raman shift is observed, which is expected for franckeite crystal due to the incommensurate stacking of the Q-layer and H-layer 22 . However, the intensity variations in Raman modes are observed as the flake thickness changes, which can be attributed to the optical interference and the associated optical field enhancement within the layered materials 22 . Next, the angle-resolved polarized Raman spectroscopy is performed to identify the in-plane rippling direction along the b-axis of franckeite crystal. Figure 2a,b show the Raman intensity contour maps depending on the linear polarization angle of the excitation beam for the 138 nm-thick franckeite flake in the parallel and perpendicular configurations, which clearly shows that the Raman mode intensities are highly polarization sensitive. The periodic variations of Raman mode intensities with respect to the linear polarization angle of the excitation beam are observed. In general, the Raman mode intensity can be expressed as where e i and e s denote the unit polarization vectors for the incident and scattered light, and R is the Raman tensor for the probed crystal. In accordance with the experimental configuration, e i and e s are in the x-y plane and there is no contribution from the z-axis. Hence, the unit polarization vectors can be written as e i = (cosθ, sinθ, 0) , while e s = (cosθ, sinθ, 0) and (−sinθ, cosθ, 0) for the parallel and perpendicular polarization configurations, respectively. The Raman tensor of A g modes for the triclinic franckeite crystal can be written as www.nature.com/scientificreports/ where u, v and w are the Raman tensor components. Thus, the Raman intensity of the A g modes in the parallel and perpendicular configurations can be written as follows 40 where // and ⊥ denote the parallel and perpendicular configurations. Figure 2c-  www.nature.com/scientificreports/ observations, it is confirmed that x-axis (0°) and y-axis (90°) are identified as the rippling direction (along the b-axis) and its perpendicular direction (along the a-axis), respectively.
Anisotropic THG response in franckeite thin flakes and determination of χ (3) . It has been recently demonstrated that franckeite crystal exhibits strong linear dichroism in visible light absorption due to its structural deformation 24 , but the anisotropic nonlinear optical response and the corresponding third-order nonlinear susceptibility tensor in these heterostructures are still unknown and need to be investigated. Figure 3a shows the transmission microscope image with the green THG emission from the 21 nm-thick franckeite flake excited by a 1560 nm pump laser with the spot size of 1.5 µm. The THG emission spectrum is plotted in Fig. 3b, showing the peak emission wavelength at 520 nm which is exactly one-third of excitation wavelength. To further confirm the THG emission, Fig. 3c shows the log-scale plot of the THG emission power as a function of the pump power, where the cubic power law fitting endorses the THG emission process. Next, the in-plane anisotropic THG emission from franckeite crystal is investigated by measuring the THG intensity dependence on the input linear polarization of the pump beam. The desired incidence linear polarization angle of the pump beam with respect to the rippling direction of the crystal is obtained by using a linear www.nature.com/scientificreports/ polarizer and a rotating half-wave plate placed in the excitation path before the sample. Also, a linear polarization analyzer is fixed at either 0° or 90° after the sample for measuring the x and y components of THG emission. Figure 3d-i plot the angular dependence of the THG emission power on the incident linear polarization angle with respect to the x-axis (rippling direction) for six franckeite flakes with the thicknesses of 21, 42, 97, 118, 138 and 192 nm, showing highly anisotropic four-lobe THG emission patterns. The black triangles represent the measured total THG emission power, while the red squares and blue circles correspond to the x and y components of the THG power. The maximum THG power is observed when the incident linear polarization angle is 0° along the x-axis (rippling direction) and the second maximum THG power occurs at 90° along the y-axis. The observed anisotropic THG emission from franckeite crystal originates from the symmetry breaking and generation of periodic ripples in franckeite lattice, which is induced by the lattice incommensurability and in-plane spatial modulation of vdW interaction between the consecutively stacked Q-and H-layers. The origin of the spatially modulated strain and generation of ripples in franckeite lattice has been recently presented 24 . The vdW adhesion energy between the consecutive Q-and H-layers is strongly dependent on the local atomic arrangements, which in turn depends on the inhomogeneous elastic deformation in the crystal lattice. As the crystal lattice minimizes its total energy (elastic plus adhesive), strong in-plane modulation of van der Waals forces and the out of plane ripples are generated in the crystal structure. As a result, the local density of states for light-matter interaction vary along both the crystal axes, leading to the anisotropic THG emission. To obtain further insights into the anisotropic THG response, the experimental observations are corroborated with the third-order nonlinear susceptibility model. As per the experimental configuration, the franckeite flake is illuminated with the linearly polarized pump beam with fundamental frequency ω. Therefore, the incident electric field can be expressed as − → E = |E| p , where p can be further connected with the rippling direction as p = xcosθ + ysinθ with θ as the linear polarization angle relative to the rippling direction (x-axis). By considering the triclinic crystal structure of franckeite, the contracted form of its third-order nonlinear susceptibility tensor χ (3) can be written as 44,45 where the first term in subscript 1, 2, and 3 denotes x, y and z respectively and the second subscript refers to the combination of three components as As the polarization of the excitation field always remains in the x-y plane in the experimental configuration, no contribution from the components of χ (3) containing z terms will be observed in the THG emission. Thus, only eight non-zero elements ( χ 11 , χ 12 , χ 18 , χ 19 , χ 21 , χ 22 , χ 28 and χ 29 ) of χ (3) will contribute to the THG emission process. Then the electric field components of THG emission and the THG intensity can be expressed as follows The obtained THG intensity expressions along x-and y-axis are further used to fit the experimentally measured values. The theoretical fitting curves are plotted as the solid curves in the respective colors in Fig. 3d-i, showing a good agreement with the measurements.
Moreover, the relative magnitudes of χ (3) tensor elements for franckeite crystal can be retrieved based on Eqs. (8) and (9).   Figure 4b plots the dependence of the THG emission power on the franckeite flake thickness for the incident linear polarization along the x-axis. The average pump power is kept at 1.1 mW corresponding to 8.65 GW/cm 2 peak irradiance. It is found that the THG emission power gradually increases up to around 4 pW at the flake thickness of 97 nm with the conversion efficiency of 3.64 × 10 -9 , afterwards it exponentially decays. Such behavior of the thickness-dependent THG emission power can be explained in terms of the optical absorption dominance on the THG emission process as the flake thickness increases. Since franckeite is a semiconductor with a narrow band gap around 0.65 eV, there is strong optical absorption of franckeite crystal at the THG emission wavelength of 3 = 520 nm (2.38 eV). The strong optical absorption of franckeite especially at large flake thickness attenuates the generated THG signal to form the observed thickness-dependent exponential decay of THG emission power. Such thickness-dependent THG signal attenuation further enables us to extract the imaginary part of the refractive index (k 3 ) at 3 = 520 nm for franckeite crystal by using the xyz 0 χ 11 cos 3 θ + χ 12 sin 3 θ + 3χ 18 cosθsin 2 θ + 3χ 19 sinθcos 2 θ χ 21 cos 3 θ + χ 22 sin 3 θ + 3χ 28 cosθsin 2 θ + 3χ 29 sinθcos 2 θ 0 x ∝ χ 11 cos 3 θ + χ 12 sin 3 θ + 3χ 18 cosθsin 2 θ + 3χ 19 sinθcos 2 θ 2 (9) I (3ω) y ∝ χ 21 cos 3 θ + χ 22 sin 3 θ + 3χ 28 cosθsin 2 θ + 3χ 29 sinθcos 2 θ , where A is a constant, l is the flake thickness. The exponential fitting of the measured THG emission power is plotted in Fig. 4b which gives k 3 = 0.83. This value agrees with the previously reported k 3 value for franckeite crystal at the wavelength of 520 nm, while the real part of the refractive index n 3 is reported as around 3.2 30 . It is noted that there is some deviation between the measured THG power and the fitted data in Fig. 4b, especially the fitted values slightly overestimate the measured ones for franckeite thin flakes with thickness less than 70 nm. One potential reason for this inconsistency is that the mechanical exfoliation of franckeite thin flakes and the transfer process from tape to substrate might influence some top surface layers, resulting in a reduced THG conversion efficiency. Another possible explanation is the surface oxidation and contamination of franckeite thin flakes exposed to air, giving a lower THG emission power. In the future work, clean mechanical exfoliation in inert gas environment and surface encapsulation to protect the pristine franckeite crystals from the surface defects, oxidation and contamination can be utilized to solve such inconsistency issue.
Taken into account the values of n 3 , k 3 and other experimental parameters of average pump power, laser pulse width, repetition rate, and spot size ( P (ω) = 1.1 mW, τ = 90 fs, f rep = 80 MHz, and W = 1.5 µm) at the pump wavelength of 1560 nm, the magnitude of χ (3) can be retrieved using the following equation 41 , where n 1 is the real part of the refractive index of franckeite crystal at the fundamental wavelength 1 = 1560 nm and here n 1 = 4 is selected, and �k = 6π 1 (n 1 − n 3 ) is the phase mismatch between the fundamental beam and the forward propagating third-harmonic emission beam in the transmission microscope optical setup arrangement. The extracted magnitude of χ (3) for franckeite crystal is 1.87 × 10 -19 m 2 /V 2 . It is worth noting that several experimental parameters in Eq. (10), especially average pump power P (ω) , spot size W, and flake thickness l due to surface roughness, may affect the measured THG signal and propagate the error in determining the χ (3) value. By considering the potential errors in average pump power P (ω) = 1.1 ± 0.1 mW, spot size W = 1.5 ± 0.06 µm and flake thickness l ± 2 nm, it is found that the extracted χ (3) value has the variation of less than ± 8%. Furthermore, while estimating the χ (3) value from Eq. (10), the real part of the refractive index at 1560 nm is selected as n 1 = 4, which is close to the reported value of franckeite in the near-infrared wavelength range with n = 4.1 at 777 nm 30 . By taking into account the uncertainty of the refractive index value selection, the n 1 value is varied from 3 to 6 in the calculation which results in the χ (3) value from 1.17 × 10 -19 m 2 /V 2 to 5.28 × 10 -19 m 2 /V 2 . It shows that the variation of n 1 does not induce an abrupt change in the χ (3) value which remains in the same order of magnitude. The extracted χ (3) value can be further justified by the previously reported χ (3) value of layered franckeite as 5.1 × 10 -11 esu (7.13 × 10 -19 m 2 /V 2 ) measured by the spatial self-phase modulation method with an ultrafast laser 23 , which is within the same order of magnitude. Additionally, in order to validate natural vdW heterostructures as an attractive alternative for advanced photonic and optoelectronic applications, the extracted χ (3) value for natural vdW heterostructures is compared with the previously reported values for other nonlinear 2D material flakes such as black phosphorous (1.4 × 10 -19 m 2 /V 2 ), graphene (1.5 × 10 -19 m 2 /V 2 ), and MoS 2 (2.9 × 10 -19 m 2 / V 2 ) 41,46,47 , showing that the χ (3) value for franckeite crystal has the similar magnitude as other 2D materials.

Discussion
In summary, we have demonstrated how the angle-resolved polarized Raman spectroscopy and the polarizationdependent THG emission can be utilized in probing the anisotropic natural vdW heterostructures with crystal structural deformation and identifying the direction of symmetry-broken periodic ripples induced by the lattice incommensurability. Both the highly anisotropic Raman response and THG emission pattern are observed from mechanically exfoliated franckeite thin flakes. We have discussed how the lattice incommensurability and inplane spatial modulation of vdW forces between the consecutively stacked PbS-like and SnS 2 -like layers play a significant role in determining the anisotropic nature of nonlinear optical process in franckeite crystal. In order to further substantiate the experimental observations, the measured THG response is corroborated with the theoretical nonlinear susceptibility model to determine the relative magnitudes of the χ (3) tensor elements and the value of χ (3) for franckeite crystal. We have also discussed how the THG emission power from franckeite flakes is a cumulative effect of both nonlinear THG emission and linear optical absorption processes, which can be utilized for tailoring nonlinear light-matter interactions in multifunctional photonic integrated circuits. The results discussed here not only provide the comprehensive understanding of nonlinear optical processes in anisotropic natural vdW heterostructures, but also pave the way to the realization of integrated quantum photonics devices and circuits. Additionally, the electrically tunable THG signal modulation has been demonstrated 48,49 , where the order of magnitude of THG conversion efficiency in graphene is tuned by controlling the Fermi energy and the incident photon energy. Beyond graphene and transition metal dichalcogenides like MoS 2 and WSe 2 , the emergence of anisotropic vdW materials has facilitated an additional degree of freedom to tailor the nonlinear optical interactions based on intrinsic structural anisotropy and crystal symmetry, which will be further harnessed to realize future design of vdW heterostructure-based anisotropic nonlinear optoelectronic devices, such as the electrically tunable broadband frequency converters for applications in optical communication, secured data transmission, encrypted signal processing, and quantum information processing.

Methods
Sample preparation. The glass substrates with 1 cm × 1 cm are treated with deionized water and isopropanol followed by ultra-sonication for 10 min to remove the undesirable residues from the surface. The processed glass substrates are used for transferring the franckeite thin flakes, which are mechanically exfoliated using Nitto tape (SPV 224) from bulk natural franckeite mineral (from San José mine, Oruro, Bolivia).
Optical setup. For collecting the Raman spectra, the sample is excited with a 632.8 nm He-Ne laser source using a 40× objective lens (NA = 0.65) and the back-reflected light is collected with the same objective lens. The incident linear polarization of laser beam is controlled using a linear polarizer and a rotating half-wave plate in the excitation path. The collected light is routed to a spectrometer (Horiba, iHR 520) using a beam splitter. An appropriate edge filter (Semrock, LP02-633RE-25) is engaged in the collection path to reject the laser light. The collected signal is passed through a linear polarizer analyzer in the collection path to further resolve the parallel and perpendicular components of Raman spectra. For collecting the THG signal, the sample is pumped with a femtosecond laser source at the wavelength of 1560 nm (pulse width 90 fs, repetition rate 80 MHz) using a 40× objective lens (NA = 0.65). The transmission light is collected through a 100× objective lens (NA = 0.7). The pump laser light is filtered out by introducing a shortpass filter in the collection path. The transmitted THG signal is then routed towards a spectrometer and an imaging charge-coupled device (CCD) camera for the measurements. www.nature.com/scientificreports/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.