Anisotropic optical responses of layered thallium arsenic sulfosalt gillulyite

Multi-element two-dimensional (2D) materials hold great promise in the context of tailoring the physical and chemical properties of the materials via stoichiometric engineering. However, the rational and controllable synthesis of complex 2D materials remains a challenge. Herein, we demonstrate the preparation of large-area thin quaternary 2D material flakes via mechanical exfoliation from a naturally occurring bulk crystal named gillulyite. Furthermore, the anisotropic linear and nonlinear optical properties including anisotropic Raman scattering, linear dichroism, and anisotropic third-harmonic generation (THG) of the exfoliated gillulyite flakes are investigated. The observed highly anisotropic optical properties originate from the reduced in-plane crystal symmetry. Additionally, the third-order nonlinear susceptibility of gillulyite crystal is retrieved from the measured thickness-dependent THG emission. We anticipate that the demonstrated strong anisotropic linear and nonlinear optical responses of gillulyite crystal will facilitate the better understanding of light-matter interaction in quaternary 2D materials and its implications in technological innovations such as photodetectors, frequency modulators, nonlinear optical signal processors, and solar cell applications.

www.nature.com/scientificreports/ ranging from energy storage to security purposes such as thermoelectrics, laser frequency modulators, radiation sensors, and photovoltaics [23][24][25] . Therefore, the insightful understanding of light-matter interaction associated with the structural properties of these materials will be valuable for their engagement in technological improvisation.
To address these issues, here we explore the anisotropic linear and nonlinear optical responses of mechanically exfoliated gillulyite thin flakes. It is noted that although CVD-grown thin films have been developed to a level of centimeter-scale size, the CVD-based growth of large-size thin films of multi-element layered materials with rational and controllable chemical compositions is still an unaddressed issue. As an alternate, naturally occurring multi-element vdW layered minerals pave an interesting way to prepare ultrathin flakes via mechanical exfoliation. However, in general the mechanical exfoliation of these naturally occurring bulk crystals in the large-size range of tens of micrometers is difficult. Here, we demonstrate that it is relatively easy to mechanically exfoliate the naturally occurring bulk gillulyite crystal into large-area thin flakes with sizes of tens of micrometers down to deep-subwavelength scale thicknesses. We further characterize the gillulyite crystal using high-resolution transmission electron microscope (HRTEM) and energy dispersive X-ray spectroscopy (EDXS) techniques to determine the structural and chemical composition information. Subsequently, the anisotropic vibrational and linear optical properties due to the low in-plane lattice symmetry are probed using polarization-resolved Raman and optical absorption spectroscopy. Furthermore, the anisotropic THG response of gillulyite crystal with respect to the incident linear polarization of pump beam is investigated and the third-order nonlinear susceptibility is retrieved from the thickness-dependent THG emission. To obtain further intuitive understanding, the thirdorder nonlinear susceptibility tensor elements are also estimated by corroborating the measured THG signal with a theoretical nonlinear model. These results will facilitate the fundamental understanding of light-matter interaction in quaternary 2D materials and advance future technological innovations in photonics and optoelectronics such as photodetection, frequency conversion, nonlinear optical signal processing, and photonic circuits.

Results
Structural and chemical composition determination of gillulyite crystal. Gillulyite belongs to the monoclinic crystal system with space group P2/n and the unit cell dimensions of a = 9.584 Å, b = 5.679 Å, c = 21.501 Å, α = γ = 90° and β = 100.07°. Figure 1a,b present the side views of the crystal structure of gillulyite projected along the b-axis and a-axis. The crystal structure of gillulyite can be described by two alternating Tlcontaining A layer and PbS-like B layer along the c-direction [26][27][28] . In the Tl-bearing A layer with partly zeolitic character, TlS 5 and As 2 S 5 groups alternate regularly along the b-direction. The B layer can be described as a substantially distorted PbS-like motif which is periodically twinned by insertion of another (As,Sb) polyhedron, and is filled with lone electron pairs and weak interactions. The order-disorder phenomena happen along the b-direction which is caused by ambiguity in the position of TlS 5 -As 2 S 5 sequences. Complex loop-branched zigzag (As,Sb) 4 S 7 double-chains of AsS 3 pyramids are formed along the a-direction in the B layer of gillulyite. It is noted that the crystal structure of gillulyite is homeotypic with those of gerstleyite Na 2 (Sb,As) 8 S 13 ⋅2H 2 O 29 and ambrinoite (K,NH 4 ) 2 (As,Sb) 8 S 13 ·H 2 O 30 . The weak interactions in the B layer of gillulyite will lead to a small interlayer cohesive energy, indicating the suitability for mechanical exfoliation of the bulk mineral. Figure 1c shows a picture of a gillulyite mineral rock, where an aggregate of vitreous deep-red platy crystals of gillulyite with a large stepped cleaved surface is presented on massive white barite. Further, the magnified view of deep-red platy crystals of gillulyite is displayed in Fig. 1d. We create the gillulyite thin flakes with different thicknesses on glass substrate via the mechanical exfoliation of the bulk natural gillulyite mineral (from Lulu Cut, South Mercur Pit, Mercur District, Oquirrh Mountains, Tooele County, Utah, USA) with Nitto tape (SPV 224). Figure 1e,f show the atomic force microscopy (AFM) images of two prepared gillulyite flakes with the thicknesses of 23 nm and 34 nm, which validate that gillulyite bulk crystal is exfoliable down to few-layer thick flakes. Large-area thin flakes with high aspect ratios can be obtained by cautiously performing the mechanical exfoliation process. Figure 1g,h display the optical reflection and transmission microscope images of two largearea thin flakes, as highlighted in region 1 with lateral dimension of 12 µm and thickness of 56 nm, as well as in region 2 with lateral dimension of 14 µm and thickness of 77 nm, where the flake thickness and surface smoothness are confirmed by the captured AFM images shown in Fig. 1i,j. Figure 1k-n display two more large-area thin flakes with thicknesses of 89 nm and 99 nm, which reaffirm the consistency of higher aspect ratios attained in the exfoliated gillulyite flakes.
Next, the gillulyite crystal is characterized using HRTEM and EDXS techniques to investigate the structural and chemical composition information. Figure 2a shows a HRTEM image of gillulyite crystal with the determined lattice spacings of ~ 4.72 Å and ~ 2.84 Å and the intersection angle of ~ 97°, which are consistent with the [200] and [020] sets of planes for the monoclinic crystal. The selected area electron diffraction (SAED) pattern in Fig. 2b further confirms the crystalline nature of the exfoliated gillulyite flakes where the spot patterns from the surface normal to the [001] crystal zone axis are displayed. We further quantify the chemical composition of gillulyite crystal. The recorded average EDXS spectrum is shown in Fig. 2c, whereas Fig. 2d-h shows the dark-field TEM image of the scanned region and the corresponding TEM-EDXS elemental maps, emphasizing the homogeneous distribution of thallium (Tl), arsenic (As), antimony (Sb) and sulfur (S) as the prime components in gillulyite crystal. In addition, the presence of copper (Cu) peak in EDX spectrum can be attributed to the TEM grid. The obtained elemental composition (atomic weight %) is summarized in Table 1, which is used to compute the compositional stoichiometry analysis, resulting an empirical formula of Tl 7.89 As 7.80 Sb 1.00 S 13.00 . The accuracy of EDXS estimation method depends upon several important factors such as the nature of probed specimen, the overlap in X-ray emission peaks, and the detection efficiency. Noticeably, the overlap in X-ray emission peaks is a crucial factor under the same operating condition. In our measurements, the X-ray emission peak of thallium is observed at 2.27 keV, whereas the corresponding sulfur peak is detected at 2.31 keV. In the EDXS spectrum, the elemental yield is estimated as per area under the curve for the respective emission peak.
Since these two peaks are significantly overlapping with each other, the estimated quantifications of thallium and sulfur with EDXS may not be precise. We anticipate that this is the major reason for the overestimation of thallium abundance and the underestimation of sulfur in the probed specimen. As a result, the empirical formula determined from the EDXS spectrum is not the excellent match with the generic gillulyite chemical formula of Tl 2 (As,Sb) 8 S 13 22 . Besides, small deviation between these two formulas can be attributed to their geological origins and the presence of carbon (C) and oxygen (O) impurities in naturally occurring minerals. Moreover, www.nature.com/scientificreports/ such deviation between the empirical and generic formula have been observed in previous studies on naturally occurring 2D materials as well 31,32 .
Polarization-resolved Raman spectroscopy of gillulyite flake. In order to understand the vibrational modes in gillulyite flakes, the polarization-resolved Raman spectroscopy is conducted on the 89 nm-thick flake in parallel polarization configuration. Figure 3a shows the recorded Raman spectra with the flake excited using a 632.8 nm He-Ne laser, while the analyzer in the collection path is set to be in the parallel direction with respect to the input linear polarization. The Raman spectrum of gillulyite crystal shows a series of distinct Raman modes within the 50-420 cm -1 frequency range, which has evident spectral similarities to the Raman spectra from other thallium arsenic sulfosalts such as lorándite TlAsS 2 , fangite Tl 3 AsS 4 and rebulite Tl 5 Sb 5 As 8 S 22 , as well as close to the Raman spectrum of orpiment As 2 S 3 , arising from the presence of AsS 3 trigonal pyramids in the crystal structures 33,34 . The Raman bands between 400 and 220 cm -1 include the characteristic four vibrational modes from ν 1 to ν 4 of the isolated and interconnected AsS 3 pyramids with band energies of ν 1 > ν 3 > ν 2 > ν 4 . The modes in the 400-330 cm -1 region represent the symmetric and antisymmetric As-S stretching vibrations of ν 1 and ν 3 , which might overlap to each other or are very close in energy. The modes in the 325-225 cm -1 region are the S-As-S bending vibrations of ν 2 and ν 4 , where the out-of-plane bending vibrations occur above 290 cm -1 and the in-plane bending vibrations appear in the low energy range. The Raman bands below 220 cm -1 originate from the lattice vibrational modes.
To obtain further insight of anisotropic vibrational modes, the dependence of the Raman modes on the incident linear polarization angle is systematically studied. Gillulyite crystal belongs to the monoclinic crystal family. Taken into considering the crystal symmetry, the Raman tensor for different A g modes can be written as 35  www.nature.com/scientificreports/ www.nature.com/scientificreports/ where a, b, c, and d are the amplitudes of Raman tensor elements, and ∅ a , ∅ b , ∅ c , and ∅ d , denote the associated phases. In addition, the Raman modes also show the significant dependence on the unit polarization vectors of the incident and scattered light e i and e s , with the Raman intensity of I ∝ |e i .R.e s | 2 . According to the experimental configuration, e i = (cosθ , sinθ, 0) and e s = (cosθ , sinθ, 0) for parallel polarization configuration, giving the Raman intensity for A g modes as follows 35 where // denotes the parallel polarization configuration and ∅ ba = ∅ b − ∅ a denotes the associated phase difference. Figure 3b shows the color map of Raman intensity as per the incident linear polarization angle for the 89 nm-thick flake. It is observed that the Raman modes are anisotropic in nature and the Raman intensity vary periodically as per the incident polarization angle. To obtain further insight on the periodicity of the Raman modes, the experimentally recorded polarization-resolved Raman spectra are theoretically fitted using Eq. (2). Figure 3c-h show the polar plots of the Raman modes at 63, 138, 187, 340, 359 and 394 cm -1 . The experimental data are denoted as black squares, whereas the theoretical fitting is shown with red solid lines, indicating a good agreement with each other. All the A g modes show anisotropic two-lobe patterns with the Raman intensity maxima either at 0° and 180° or 90° and 270°. Moreover, the orientation of A g modes in parallel polarization configuration indicates the crystal axis of the investigated gillulyite flake 13,36 . Here, the 0° and 90° directions are acknowledged as the a-axis (along x-axis) and b-axis (along y-axis) of gillulyite crystal, for further analyzing the anisotropic linear and nonlinear optical properties of the crystal. It is worth noting that the effects of glass substrate on the observed anisotropic Raman modes from gillulyite crystals are negligible (see Supplementary  Information). Furthermore, the Raman spectra of four different gillulyite crystals with thicknesses of 23, 56, 77, and 89 nm are collected under the excitation of linearly-polarized 632.8 nm laser source without the analyzer in the collection path. The collected Raman spectra are shown in Fig. 4a and the Raman A g peaks at 63, 138, 187, 340, 359, and 394 cm -1 are indicated with black dashed lines. It is observed that the Raman signal intensity increases with the increased gillulyite crystal thickness. However, there is no significant frequency variation in the observed Raman modes as a function of the crystal thickness. Next, we perform the polarization-dependent Raman spectroscopy of these crystals in parallel polarization configuration and compare the anisotropy ratios of the Raman intensities for the observed A g modes at 63, 138, 187, 340, 359, and 394 cm -1 . The anisotropy ratio depending on the crystal thickness is plotted in Fig. 4b. It is observed that most of the Raman A g modes are highly anisotropic in nature and the anisotropy ratio just slightly changes as per the probed gillulyite crystal thickness.
Polarization-resolved optical absorption spectroscopy of gillulyite flake. The prevalent anisotropy in polarization-dependent Raman modes also accredit the presence of linear dichroism in gillulyite crystal. To get further intuition regarding the intrinsic linear dichroism resulted from the low in-plane crystal symmetry, the optical absorption characteristics of two gillulyite flakes with thicknesses of 56 nm and 89 nm are investigated in the visible range from 450 to 800 nm using polarization-resolved optical absorption spectroscopy. Figure 5a,e,i show the measured reflectance (R), transmittance (T) and absorbance (A = 1 − R − T) spectra for these two gillulyite flakes, by keeping the incident linear polarization fixed along x-axis of the probed crystal. It is found that the reflectance is higher for the relatively thick flake (89 nm), whereas the transmittance is higher for the thinner flake (56 nm). Interestingly, the reflectance spectra from both flakes show a gradual decrement with a dip around 725 nm, whereas the transmittance increases monotonously within the range of 450 nm to 800 nm. www.nature.com/scientificreports/ The reflectance dip around 725 nm results in an absorbance peak shown in Fig. 5i, which is related to the optical band gap of gillulyite crystal, giving the mineral its deep red color similar to that of lorándite and fangite 37 .
To further investigate the anisotropic features of optical absorption and linear dichroism, polarization-resolved reflectance, transmittance, and absorbance spectra are systematically measured for the 89 nm-thick gillulyite flake, as shown in Fig. 5b,f,j. The linear polarization angle is defined relative to x-axis of the crystal. The polar plots of reflectance, transmittance, and absorbance spectra as a function of the linear polarization angle at two different wavelengths of 500 nm and 700 nm are further displayed in Fig. 5c,d,g,h,k,l. The measured data points are theoretically fitted using the formula α(θ) = α x cos 2 (θ) + α y sin 2 (θ) , where α x and α y are the magnitudes along x-axis and y-axis. The obtained reflectance, transmittance, and absorbance anisotropy ratios at 500 nm (700 nm) are R x /R y = 1.25 (1.35), T x /T y = 0.82 (0.89), and A x /A y = 0.65 (0.86), respectively. It is indicated that the gillulyite crystal absorbs photons anisotropically along the directions of a-axis (x-axis) and b-axis (y-axis), exhibiting strong wavelength-dependent linear dichroism effects. It is noted that the effects of glass substrate on the reported transmittance, reflectance, absorbance values and the linear dichroism from gillulyite crystals are negligible (see Supplementary Information).
Anisotropic third-harmonic generation response of gillulyite flake. The reduced in-plane symmetry in gillulyite crystal also suggest the high anisotropic nonlinear optical response. The anisotropic THG emission in gillulyite flakes is probed using a 1560 nm pulsed laser with the beam waist of 1.5 µm. Figure 6a shows the measured THG spectrum from the 89 nm-thick gillulyite flake with the emission peaks at 520 nm, which is one-third of the fundamental wavelength. The THG emission process is further confirmed by the cubic power law fit between the THG emission power and the incident pump power, as shown in Fig. 6b. Next, the inplane anisotropic THG emission as per the incident pump beam polarization is characterized. The desired input polarization is obtained by using a linear polarizer and a rotating half-wave plate in the illumination path. The x-and y-components of THG emission power are filtered out by introducing an analyzer orientated in parallel (0°) and perpendicular (90°) to the crystal's a-axis. The obtained THG emission patterns are highly anisotropic with four-lobe patterns. The maximum THG emission power is collected at 0° and 180° along the crystal's a-axis (x-axis), while the secondary maxima is recorded at 90° and 270° along the crystal's b-axis (y-axis). The measured x-component, y-component and total THG power data points are represented as red squares, blue dots and black upper triangles, respectively, which are further corroborated with the theoretical third-order nonlinear susceptibility model. The third-order nonlinear susceptibility tensor for a monoclinic crystal can be written as follows 38,39  www.nature.com/scientificreports/ 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 Taken the experimental configuration into account, the incident linearly polarized pump beam can be expressed as � E =x(|E|cosθ ) +ŷ(|E|sinθ) , where θ is the linear polarization angle relative to the crystal's a-axis. Also, the incident pump beam resides only in the x-y plane, therefore the third-order nonlinear susceptibility tensor elements containing z-components can be neglected. Hence, the THG intensity components can then be written as 36,40,41 Equation (4) and (5) are used to fit the measured angular dependence of THG emission. The theoretical fittings are shown as the solid curves in the respective colors in Fig. 6c-f, showing a good agreement with the measured data.
Additionally, the theoretical fittings enable us to retrieve the relative magnitudes of nonlinear susceptibility tensor elements χ 11 ,χ 18 ,χ 22 and χ 29 , as shown in Fig. 7a for different flake thicknesses. No significant variation in the values of χ (3) elements is observed for flakes with different thicknesses, and the average values are χ 11 : χ 18 : χ 22 : χ 29 = 1 : 0.119 : 0.620 : 0.131 . Furthermore, the average ratio of |χ 11 | 2 /|χ 22 | 2 indicates the THG anisotropy ratio I (3ω) x (θ = 0)/I (3ω) y (θ = 90) in the gillulyite flakes, which almost remains as a constant of 2.6. Finally, the third-order susceptibility χ (3) value of gillulyite crystal is estimated. Figure 7b plots the measured THG emission power as a function of the flake thickness. The average power of pump beam is 1.3 mW with 10.22 GW/cm 2 peak irradiance. The recorded THG emission power gradually increases up to 0.74 pW for the 133 nmthick gillulyite flake and then decays exponentially. This thickness-dependent THG emission process can be explained in the context of two competitive phenomena of optical gain and loss. In case of the thin flakes, the collected THG emission is proportional to the square of the flake thickness, so that the THG emission power increases for the flake thickness up to around 149 nm. Afterwards, for the thick flakes greater than 149 nm, the optical absorption starts playing a dominant role in attenuating the THG signal propagation, which leads to the exponential decay. The exponentially decayed THG emission further enables us to extract the imaginary part of refractive index (k 3 ) at 3 = 520 nm by fitting the thickness-dependent THG emission power equation where A is a constant, d is the flake thickness. Figure 7b shows the measured data (black squares) and the equation fitting (red curve) with k 3 = 0.549. It is worth noting that the real parts of refractive indices of gillulyite correspond closely to those of fangite, in which n 3 = 2.81 at 589 nm and n 1 = 2.60 at 1553 nm 37,42  www.nature.com/scientificreports/ With the real parts of refractive indices of gillulyite crystal n 1 = 2.60 at 1 = 1560 nm and n 3 = 2.81 at 3 = 520 nm, the estimated magnitude of χ (3) for gillulyite crystal is 2.059 × 10 -20 m 2 /V 2 .

Discussion
In summary, we have demonstrated the mechanical exfoliation of large-area thin gillulyite flakes of various thicknesses. The structural and chemical composition information of gillulyite crystal has been characterized using HRTEM and EDXS techniques. The anisotropic linear and nonlinear optical properties of gillulyite thin flakes due to the reduced in-plane crystal symmetry have been investigated. The anisotropic Raman modes and linear dichroism have been probed using polarization-resolved Raman and optical absorption spectroscopy. The anisotropic THG emission process in gillulyite crystal has also been explored with the extracted third-order nonlinear susceptibility. In the line of these findings, our presented discussion of gillulyite crystal in the context of anisotropic structural, vibrational, and optical properties will facilitate an insightful comprehension of the light-matter interaction in thallium-contained complex 2D materials as well as pave the way for its implications in photodetection, laser frequency conversion, nonlinear optical signal processing, encrypted optical communication, and photonic circuits.
Finally, the uniqueness in gillulyite crystal structure facilitates an exceptional platform for the wide range tunability of vibrational, optical and electronic properties by comparing with its component binary sulfides via stoichiometric engineering and heteroatoms doping, which will advance many technological innovations. As one thallium-containing sulfosalt compound, gillulyite can be identified as a potential 2D material for realizing several unique applications ranging from thermoelectrics, frequency modulators, to spintronics and radiation sensors. The demonstrated mechanically exfoliated large-area thin gillulyite flakes offer a convenient platform for ultracompact multifunctional flatland optical, electronic and optoelectronic integrated devices, because rational and controllable synthesis of such complex 2D materials remains a challenge till date. Besides, the demonstrated polarization-sensitive vibrational and optical responses in gillulyite crystals indicate that the input polarization can be harnessed as an optical switch in applications such as phonon-based frequency modulation, polarized photodetection, encrypted data transfer, and nonlinear signal processing. With this hindsight, we anticipate that vdW layered sulfosalt mineral gillulyite truly has the potentials to be applied into future integrated photonic, electronic and optoelectronic technologies.

Methods
Sample preparation. The sample preparation process involves two steps. In the first step, the glass substrate with 1 cm × 1 cm is treated with deionized water, acetone and isopropanol, followed by ultra-sonication and dried with N 2 gas. This pretreatment of the glass substrate is repeated several times to ensure the removal of undesired residues on the glass substrate. Further, the glass substrate is heat treated at 150 -170 °C for 5 min to remove any remnant solvent on the glass surface. In the second step, gillulyite flakes are mechanically exfoliated using Nitto tape (SPV 224) from naturally occurring bulk gillulyite mineral (from Lulu Cut, South Mercur Pit, Mercur District, Oquirrh Mountains, Tooele County, Utah, USA). As the mechanical exfoliation process is a top-down method, it requires multiple efforts to attain the flake thickness of sub-hundred nanometer scale. After completing the exfoliation process, the gillulyite thin flakes are transferred to the pretreated glass substrate followed by heat treatment at 120-130 °C for 2 min. Then, the exfoliated flakes are examined using optical reflection and transmission microscope, and atomic force microscope for estimating the flake thickness. This procedure is repeated multiple times till the desired flake thickness range is obtained. In the view of this, the mechanical exfoliation process can certainly control the thickness range of the prepared flake.
Polarization-resolved Raman spectroscopy. The sample is illuminated with a 632.8 nm He-Ne laser using a 40 × objective lens (NA = 0.65) and the back-reflected signal is collected to a spectrometer (Horiba, iHR 520) with a beam splitter. The polarization of incident beam is controlled by introducing a linear polarizer and a rotating half-wave plate in the excitation path. The elastic scattered light is filtered out by engaging an edge filter (Semrock, LP02-633RE-25) in the collection path. The collected signal is further passed through a linear polarization analyzer to record the parallel polarization component of the Raman signal.
Polarization-resolved absorption spectroscopy. A broadband white light source (Thorlabs, SLS201L, 360-2600 nm) is passed through a linear polarizer and focused on the sample with an 80 × objective lens (NA = 0.5). To record the reflection spectra, the back-reflected light is routed to the spectrometer using a beam splitter, whereas in the case of the transmission spectra, the transmitted light through the sample is collected using another 100 × objective lens (NA = 0.7). The reflection and transmission spectra are normalized with respect to the light source spectrum to obtain the reflectance (R) and transmittance (T) spectra, so that the absorbance (A) spectrum is obtained with A = 1 − R − T.
Third-harmonic generation measurement. The sample is pumped with a femtosecond laser source at the fundamental wavelength of 1560 nm (repetition rate 80 MHz, pulse width 90 fs) using a 40 × objective lens (NA = 0.65). The transmitted THG signal through the sample is collected by a 100 × objective lens (NA = 0.7). The  www.nature.com/scientificreports/ pump beam is rejected by introducing a shortpass filter in the collection path. The THG signal is then routed to a spectrometer (Horiba, iHR 520) and an imaging camera for recording the spectra and corresponding image.