Spin-selective orbital reconfiguration and colossal nonlinear anisotropy in defect-engineered atomically thin quantum dots

Abstract

Spin-texture and wavefunction modification is still a central problem in the attempt to control spin orientation and distribution, realize effective spin transport, minimize energy dissipation, and improve functionality in spin-based information processing. Here, we studied a specific combination of atomically thin, nonmagnetic, group-13-based post-transition metal chalcogenide (PTMC) quantum dots (M₁₀X₁₂; M = Ga, In, Tl; X = S, Se, Te) with a single transition metal (Zr, Mo, Mn) introduced at the central site using density functional theory. It is observed that the spatial distributions of frontier orbitals are not only spin-selective but also site-selective in real-space due to the induced TM–X₃ hybridization, which could allow independent tuning of both spin and spatial characteristics simultaneously. The normalized density overlap can differ by up to 10⁻³, whereas the normalized signed-amplitude overlap of the corresponding wavefunctions can reach values less than or equal to 10⁻⁴. The local and effective magnetic moments of the M₉TMX₁₂ structure can range from 1 → 7µ_B and 1.62 → 20.93µ_B, respectively, while the spatial extent of identical-spin orbitals can differ by up to ∼197 a.u. Moreover, the energy gap of the pristine combination ranges between 0.89 and 5.02 eV, and it widens upon TM substitution for E^↑_g (2.6–5.41 eV) while narrowing for E^↓_g (1.62–5.36 eV), which demonstrates that one spin channel can consistently remain energetically more accessible. A clear directional imbalance in the nonlinear optoelectronic response is observed, with the induced polarization switching between in-phase and out-of-phase (phase inversion) based on the field direction. The hyperpolarizability components reach maximum values of up to ∼10⁶ (β_‖) and ∼10⁹ (γ_‖, γ_⊥), while the relative distribution between parallel and perpendicular directions remains consistent. Such multifunctional responses resulting from a single-site TM perturbation facilitate fundamental insight into how localized electronic changes simultaneously modulate spin selectivity, orbital anisotropy, and nonlinear polarization as interdependent quantum variables in low-dimensional nanoflakes, which holds promise for the design of multivariable quantum information processing architectures and next-generation nano-spintronic materials.

---

1 Introduction

The convergence of low-dimensional quantum materials and spintronics serves as a practical basis to investigate spin-dependent phenomena at the atomic scale. At present, the traditional semiconductor chips heavily rely on bulk materials to achieve charge quantization of carriers for binary logic, rather than spin quantization.^1,2 Unlike bulk semiconductors with repeating crystal lattices, the two-dimensional (2D) nanoflakes can offer greater efficiency and controllability of electrons due to their finite-scale, surface-dependent, and non-periodic nature.³ In contrast to the quasi-continuous energy bands exhibited by bulk crystals, 2D nanoflakes referred to as quantum dots (QDs) show discrete energy states which leads to tunable energy gaps,^4,5 sharp absorption/emission,⁶ edge-dominance⁷ and precise carrier localization.^8,9 Perhaps 2DQDs are currently among the strongest candidates for controlling spin textures, placing them as a vital platform for creating novel spintronic and quantum computational devices.^10,11 Influenced by the atomic sheets of graphene-like materials,¹² several classes of other layered materials have been evaluated for different applications. Beyond transition metal dichalcogenides (TMDs), a new set of materials named post-transition metal chalcogenides (PTMCs) have gained enormous interest due to their exceptional polymorph-dependent properties, surface chemistry, and optoelectronic and photovoltaic applications.^13–18 However, previous studies on PTMCs to investigate magnetism include doping induced ferromagnetic (FM) ordering,¹⁹ magnetization of nonmagnetic materials,²⁰ vdW magnetism,²¹ and magnetotransport²² in bulk or extended 2D structure. Moreover, defect-induced magnetic semiconductors have been studied to understand Curie temperature variation as a function of impurity concentration in PTM-based III–V semiconductors,²³ along with strain-tunable magnetic semiconductors for magnetic anisotropy²⁴ and magnetic proximity effects in non-magnetic–ferromagnetic bilayers.²⁵ Although much has been achieved, modulation of doping-induced electronic features or redistributions of wavefunction at the single atomic layer scale for much greater controllability has not yet been studied thoroughly, as we require a careful strategy and systematic approach to achieve such subtle modifications.

Previous studies on wavefunction modification in 2D or dot-based systems include experimental efforts to understand the electronic structure, charge localization, and optical properties induced with different point defects,^26–28 as well as the manipulation of orbital degrees of freedom through inter-dot coupling, dot deformation, or external fields, using microscopy and DFT calculations.^29–31 However, these studies do not answer how the spatial distributions or defect-induced relocalizations of orbitals contribute, which could have built a more constructive understanding of charge localization mechanisms. Present works^32–34 however address primarily problems and challenges such as electronic wave functions to control optical and electronic properties, orbital quenching from strong crystal field effects, wavefunction overlap in multi-dot arrays, which leads to a loss of well-defined orbitals, and lack of independent control on spin and orbital states. Recently, it has been shown that, unlike conventional TMDs, PTMC-based thin flakes have heavy elements with large atomic masses, which enhance orbital contributions to magnetism and make them a suitable platform to study orbital-dominated and anisotropic magnetic responses.³⁵ Because their d-orbitals are fully occupied and chemically inert, PTM-based low-dimensional systems tend to form s–p orbital interactions when bonded with chalcogens. Due to such characteristics, their ability to respond to point-defect perturbations would likely be more pronounced and the modulation of electronic wavefunctions would be more effective under atomic confinement. Substitutional doping, being one of the earliest and widely used perturbative approaches in semiconductor technology, has long been used to alter the physicochemical properties of materials by introducing localized perturbations in the potential landscape.³⁶ In particular, the partially filled d-orbitals of TM atoms can favor manipulation of carrier concentration, Fermi level positioning, defect-induced states, band alignment, and band offsets within the host lattice. In this context, p–d orbital hybridization, by means of substitution of group-13 cations, is expected to form localized, symmetry-broken electronic states in real-space in confined 2D metal chalcogenide nanodots. This can minimize orbital quenching and allow greater spatial and spin-orbital control under weak crystal field constraints.

Magnetoelectronic properties are not the only physicochemical characteristics modified by doping, as both linear and nonlinear optical responses are also effectively reconfigured in the nanodots.^37,38 A previous study reveals how metal substitution in 2D α-In₂Se₃ semiconductors creates Janus structures like InGaSe₃ and GaInSe₃³⁹ which results in enhancement of nonlinear optical properties, particularly second-harmonic generation (SHG), making them a promising candidate for advanced photonic devices. Furthermore, recent spectroscopic measurements found that quantification of third-harmonic generation (THG) is of great importance for layered III–VI compounds as nanoflakes of GaSe and GaTe showed strong third-harmonic peaks at wavelengths of ≈520 nm under an excitation wavelength of 1560 nm.⁴⁰ Such strong nonlinear optical responses, in both static and dynamic limits, are important for developing applications in electro-optic modulations, optical nano-sensors, and optical communications due to their higher-order susceptibility.^41–43 Hence, to build a complete description of such interrelated aspects, a constructive theoretical study is essential to account for how a small defect would govern electronic disposition, spatial distribution of spins, and nonlinear polarization in a confined molecular region.

This work presents a framework to investigate how modification of only a single site in nonmagnetic PTMC (M₁₀X₁₂; M = Ga, In, Tl and X = S, Se, Te) nanoflakes can lead to the independent control of spatial localization and higher-order hyperpolarizabilities, whose sizes are comparable to the sub-de Broglie length scale. Besides wavefunction relocalization and orbital orthogonality, we found large anisotropic g-factor shifts, high spin-polarization, and spin-selective trap-like localized states, which appear to be major contributors for having enormous directional linear and nonlinear optoelectronic responses, as observed from the partial density of states (PDOS). The isosurface plots of the highest occupied and lowest unoccupied molecular orbitals (HOMO/LUMO), along with calculation of the electronic spatial extent (ESE), reveal that perturbed wavefunctions can simultaneously exhibit both localization and delocalization characteristics relative to the substituted site, mainly controlled by the unpaired d-electron count (4d², 4d⁵, 3d⁵). In addition, some QDs exhibit phase-inversion characteristics with a consistent ratio between first- and second-hyperpolarizability relative to the induced polarizations. Overall, these findings show how effectively a nonmagnetic host system can control the intrinsic spin-dependent spatial characteristics over its confined region for possible applications in spintronics, magneto-optical modulation, reconfigurable quantum logic gates, quantum sensing, and magnetic and spin memory technologies.

2 Computational method

Multifunctional dependencies of the considered nanodots were calculated within DFT approximation and implemented in the Gaussian 16W program.⁴⁴ The magneto-electronic and optoelectronic calculations were performed with a range-separated hybrid exchange–correlation functional ωB97X-D⁴⁵ and an effective core potential basis set LANL2DZ.⁴⁶ To accurately describe the spatial distributions of molecular orbitals in a non-periodic finite size structure, it is important to apply such theories that are able to include long-range exact exchange, reduce self-interaction errors, as well as improve the descriptions of spin-dependent orbital interactions.^47,48 Both closed-shell and open-shell approaches were used to perform geometry optimization on the systems with the time-independent self-consistent field (SCF) method to confirm different spin-polarized states. In accordance with Hund's rule, the ground state of the substituted nanoflakes corresponds to a spin-polarized configuration with the highest possible multiplicity, as it minimizes the total energy. All triangular QDs were fully optimized on the potential energy surface until the maximum force converged below 1.5 × 10⁻⁵ Hartree Bohr⁻¹ and the RMS change in the density matrix fell below 1.5 × 10⁻⁸ Hartree Bohr⁻¹. Spin–orbit coupling effects were captured perturbatively through the orbital-Zeeman and g-tensor contributions. The g-tensors were determined with the gauge-including atomic orbital (GIAO) method, and linear/nonlinear optoelectronic tensors were calculated with the finite-field method. The detailed computational methods are provided in the SI.

3 Results and discussion

The discussion of magneto-electronic and nonlinear optoelectronic responses is organized into five interconnected parts: (i) core electronic trends, (ii) frontier orbital characteristics, (iii) orbital spatial extent, (iv) magnetic strength and spin-dependent response, and (v) static optical response. In the first subsection, we discuss the systematic, controllable electronic features observed in the considered nanodots. In the following subsection, the rearrangement and redistribution of the molecular orbitals are discussed within a broader scientific context. In the next subsection, we describe the resultant intrinsic characteristics and patterns of the orbital wavefunctions in both pristine and centrally substituted configurations. The discussion is then continued to address the observed magnetic responses and field-induced magnetic characteristics. Finally, our discussion comes to an end by highlighting notable nonlinear components and their hybridization-dependent characteristics across different substituted configurations, followed by a concise summary of a few possible ways to validate these results experimentally.

Before going into the detailed discussion, an overview of the figures is presented briefly. Fig. 1 and 2 show the atomic structure of the quantum dots discussed in this work. Fig. 3(a)–(c) describes the projections of the total contributions from each orbital's density of states in the quantum dots. Fig. 4 represents the HOMO and LUMO distributions of the pristine configurations, while Fig. 5(a)–(i) depicts the spin-dependent HOMO–LUMO distributions in the centrally perturbed configurations. Fig. 6 shows how the spin densities are spread across different spatial regions. Finally, Fig. 7 compares the directional uniformity preserved in the systems, even after an increase in their respective anisotropic components. The following subsections describe these figures in detail, using a topic-based order that best supports the flow of the discussion rather than a strict numerical order.

Fig. 1 (a)–(i) Single-atom-thick, triangular M₁₀X₁₂ quantum dots with zigzag edges, structurally derived from honeycomb-shaped 2D metal chalcogenide nanoflakes. Each bond length ranges from approximately 2.1 Å to 3.2 Å. The edge-effects and geometric confinement have a critical influence on their electronic structure and emergence of distinctive properties.

Fig. 2 (a)–(i) TM-substituted counterparts of the pristine nanodots with chalcogen-coordinated point defects (X₃ sites). Such localized defects disrupt the inversion symmetry of the configurations and demonstrate highly controllable multifunctional responses in real-space.

Fig. 3 (a) Partial density of states plotted for different orbitals in Ga-based systems. The data correspond to pristine Ga₁₀X₁₂, Zr-doped Ga₉S₁₂, and Mo-doped Ga₉(Se,Te)₁₂ configurations, with internal titles indicating the specific composition in each case. (b) Partial density of states plotted for different orbitals in In-based systems. The data correspond to pristine In₁₀X₁₂, Zr-doped In₉S₁₂, and Mo-doped In₉(Se,Te)₁₂ configurations, with internal titles indicating the specific composition in each case. (c) Partial density of states plotted for different orbitals in Tl-based systems. The data correspond to pristine Tl₁₀X₁₂ and Mn-doped Tl₉(S,Se,Te)₁₂ configurations, with internal titles indicating the specific composition in each case. (a)–(c) Orbital-projected PDOS plotted for the pristine and the doped configurations comprises each orbital's total contribution. All the plots show a well-defined HOMO–LUMO gap with distinct lower energy and higher energy spectra in the potential landscape. For the pristine QDs, the largest HOMO–LUMO gap was derived for Tl₁₀S₁₂ (5.02 eV) and the lowest was found for Ga₁₀S₁₂ (0.89 eV). After the central perturbation with TM, emergence of new hybridized electronic states around and near the energy gap is observed. These localized states are spin-dependent by nature, as we can see that the overlapping states energetically vary between the alpha- and beta-spin orbitals. Such variation is very high when TM's d-orbital has higher unpaired electrons, and this controllable feature not only allows modification of pristine's orbital redistributions and fine texture of spin densities, but also enhances the nonlinear responses and the intrinsic spin–orbit coupling. The largest energy gap was exhibited by α-In₉ZrS₁₂ (5.41 eV) while the smallest was captured for β-Ga₉MoTe₁₂ (1.62 eV). All PDOS are aligned to a common reference energy (set at 0 eV) for each system to allow direct comparison.

Fig. 4 Isosurface plots of the frontier molecular orbitals for pristine configurations. Red and green regions represent the phase of the orbital wavefunctions, showing where electron density is distributed between atoms. This distribution tells us at which location in real-space, bonding interactions of electrons occupy (the most probable) a molecular region when just below (HOMO) or above (LUMO) the energy gap center. Such asymmetric lateral distribution of MO explains the intrinsic electronic structure of the nonmagnetic host and provides a baseline to understand how doping-induced perturbations would modulate orbital character in a site-specific manner.

Fig. 5 (a)–(i) Isosurface plots of the frontier molecular orbitals for TM-doped configurations. Orbital relocalization and distinct spin-selective distribution are observed due to the formation of local TM–X₃ point defects. The wavefunctions show spin-dependent localization, either concentrated around the defect center or redistributed toward edge sites, and remain spatially orthogonal for same-spin orbitals of different energies. Among the configurations, β of Zr–(S)₃ and Mn–(S/Te)₃ coordination exhibits comparatively more pronounced spatial modulation, charge concentration, and orbital reorganization in the host quantum dots.

Fig. 6 (a)–(i) Spin density distribution of spin-up (green) and spin-down (blue) orbitals in the defected nanodots. The dominant spin-down components show preferential spin-polarization with distinct spatial localization around the central defect or redistribution toward specific edge sites, depending on the host configuration. This answers how selective hybridization tends to favor specific spin channels that can have access to more electronic states (hence resulting in a lower energy gap in the corresponding orbitals) within a confined geometry.

Fig. 7 Directional sensitivity (R_A, R_B, R_Γ) derived from linear polarizability, first hyperpolarizability and second hyperpolarizability, respectively. This comparison basically describes the degree of directional uniformity retained in the substituted systems, even though they exhibit colossal nonlinear intensities along different components. R_A ranges from 0.97 to 1.15 for pristine systems and from 0.67 to 12.94 for centrally perturbed systems, which indicates that the anisotropic response of linear polarizability remains relatively uniform in most cases, but shows significant variation in Mn-doped Tl₉Se₁₂. In contrast, the static hyperpolarizabilities exhibit much more consistency in directional behavior, with all the configurations having R_B values ranging from 2.9 to ∼3, and R_Γ from 2.9 to 3.1. The only exception is the In₁₀Se₁₂ QD, which has a ratio of 0.34, indicating that its orthogonal direction is more responsive to the external field. Such preserved directional sensitivity with a high increment in nonlinear response suggests the substantial anisotropic behavior of these nanoflakes under an external electric field.

3.1 Core electronic trends

We begin this analysis by presenting the structural information and addressing the stability of the considered structures. We then describe the overall trends observed in several electronic descriptors, and conclude by discussing how the observed features exhibit clear tunable spin-dependent responses. The structural dimensions and local bonding of the pristine and central-site substituted nanoflakes were quantified using the maximum atomic extent, radius of gyration, and bond lengths. As presented in Table 1, D_max was calculated as the largest distance between any two atoms in the configuration, and R_g was computed from the root-mean-square distance of all atoms from the center of mass. The corrugation measures how much atoms deviate out-of-plane relative to the nanoflake's planar structure (i.e., buckling), which helps to understand the presence of any local distortions. These descriptors provide a quantitative measure of the overall size and compactness of the finite nanoflake structures. Across all systems, the pristine nanoflakes exhibited planar configurations with zero corrugation and average bond lengths ranging from 2.345 Å (Ga₁₀S₁₂) to 2.937 Å (Tl₁₀Te₁₂), while D_max values varied from 12.733 Å to 15.586 Å and R_g from 5.021 Å to 6.269 Å. Upon substitution, the nanoflakes maintained planarity, with only minor changes in D_max (e.g., 13.110 Å for Ga₉ZrS₁₂ vs. 12.733 Å for Ga₁₀S₁₂) and R_g (5.034 Å vs. 5.021 Å). Dopant–neighbor bond lengths were measured for each locally perturbed system, which showed values consistent with local structural relaxation: 2.523–2.558 Å for Zr–S, 2.517–2.527 Å for Mo–Se, 2.750 Å for Mo–Te, and 2.391–2.525 Å for Mn–chalcogen bonds. Overall, these metrics indicate that the substituted nanodots preserve the global size and compactness of their pristine counterparts while accommodating local bond rearrangements around the site of substitution.

Table 1 Structural parameters of pristine and TM-substituted PTMC QDs. D_max: maximum diameter; R_g: radius of gyration; d̄: average bond length; d_TM–X: dopant–neighbor bond length

System	Pristine QDs				Substituted QDs
	Dmax (Å)	Rg (Å)	d̄ (Å)	dmin–dmax (Å)	Dmax (Å)	Rg (Å)	d̄ (Å)	dmin–dmax (Å)	dTM–X (Å)
Ga-based QDs
Ga10S12/Ga9ZrS12	12.733	5.021	2.345	2.130–2.625	13.110	5.034	2.366	2.130–2.764	2.558, 2.523, 2.558
Ga10Se12/Ga9MoSe12	13.217	5.270	2.470	2.259–2.753	13.260	5.267	2.466	2.264–2.786	2.527, 2.517, 2.527
Ga10Te12/Ga9MoTe12	14.075	5.656	2.659	2.455–2.954	13.713	5.605	2.624	2.538–2.750	2.750, 2.750, 2.750
In-based QDs
In10S12/In9ZrS12	13.631	5.369	2.495	2.293–2.740	13.662	5.397	2.504	2.294–2.750	2.432, 2.432, 2.432
In10Se12/In9MoSe12	14.128	5.627	2.628	2.421–2.890	13.740	5.526	2.585	2.483–2.733	2.663, 2.637, 2.636
In10Te12/In9MoTe12	14.990	6.017	2.822	2.609–3.113	14.451	5.895	2.768	2.691–2.874	2.788, 2.788, 2.788
Tl-based QDs
Tl10S12/Tl9MnS12	14.292	5.646	2.625	2.406–2.891	14.298	5.517	2.587	2.391–2.920	2.391, 2.449, 2.391
Tl10Se12/Tl9MnSe12	14.748	5.888	2.748	2.530–3.022	14.585	5.756	2.709	2.495–3.029	2.525, 2.495, 2.525
Tl10Te12/Tl9MnTe12	15.586	6.269	2.937	2.715–3.236	15.950	6.174	2.912	2.670–3.269	2.671, 2.693, 2.670

---

In pristine [Ga/In/Tl]–[S/Se/Te] nanoflakes, the 4p, 5p, and 6p orbitals of the metal atoms primarily interact with chalcogen's 3p, 4p, and 5p orbitals, while the filled d-shells remain as spectators. The binding energies (E_b) have been derived to measure and quantify these interaction strengths with eqn (S1) (see the SI). As presented in Table 2(b), for all the MX systems, it is found that E_b increases when larger chalcogens are placed, due to their more diffuse orbitals and effect on charge redistribution. The highest (3.64 eV) E_b is observed for Ga₁₀S₁₂, while the lowest (2.53 eV) is derived for Tl₁₀Te₁₂ in pristine dots. It is also observed that after substitution of a single metal atom with a TM (Zr, Mo, Mn) from the central site, the interaction energy increases, as the TM's partially filled d-orbitals can hybridize more frequently with the orbitals of the host elements. The emergent d–p covalent framework strengthens the interatomic bonding which maintains the same chalcogen-dependent energy hierarchy. The smallest increment (0.2 eV) occurs in Mn-substituted Tl₉Te₁₂, whereas Mo-substituted Ga₉Te₁₂ exhibits the largest (0.27 eV). Conceptually, the binding energies shift scales with the electronegativity difference between the TM dopant and the chalcogen:

ΔE_b ∝ (χ_TM − χ_X)²

where χ_TM is the Pauling electronegativity of the TMs and χ_X represents the chalcogens. Due to long-range structural relaxation and electronic redistribution induced by the dopant, the replacement of the host atom by a single TM causes a local perturbation in both ionic size and electronegativity. Because these TMs have larger atomic radii and different valence p-orbital energies than the pristine metal atom, the local bonding environment around the substitution site undergoes strain and charge reorganization (as per the dipole moments and isosurface plots). This perturbation propagates along the finite quantum dot geometry, which leads to nonlocal relaxation effects that appear as subtle distortions at the edges (Fig. 2). However, such structural distortions do not lead to instability of the configurations, as confirmed by the formation energies, vibrational analysis, and energetic response under finite temperature perturbations.

Table 2 Electronic and magnetic properties of group-13-based MX configurations. (a) Pristine structures, (b) TM-substituted reconstructed structures, and (c) magnetic responses due to interaction with an external field

(a)
Host	Eb (eV)	Eg (eV)	μ0 (a.u.)
Ga10S12	3.64	0.89	0.0022
Ga10Se12	3.32	1.30	0.0026
Ga10Te12	2.83	3.50	0.0203
In10S12	3.47	3.90	0.0029
In10Se12	3.17	4.80	0.0013
In10Te12	2.72	4.60	0.0015
Tl10S12	3.15	5.02	0.0242
Tl10Se12	2.92	2.85	0.0169
Tl10Te12	2.53	4.40	0.0261

---

(b)
Host	Dp	Eform (eV)	Eb (eV)	E^↑g (eV)	E^↓g (eV)	μ0 (a.u.)
Ga9S12	Zr	−4.78	3.86	5.30	4.90	1.8587
Ga9Se12	Mo	−4.90	3.55	4.36	2.10	4.2452
Ga9Te12	Mo	−4.54	3.10	5.33	1.62	5.9322
In9S12	Zr	−5.17	3.70	5.41	5.36	0.0033
In9Se12	Mo	−4.74	3.40	5.23	1.70	6.8281
In9Te12	Mo	−4.88	2.94	5.50	1.65	0.0079
Tl9S12	Mn	−4.84	3.40	4.85	1.94	3.0787
Tl9Se12	Mn	−4.20	3.10	4.66	1.88	1.9514
Tl9Te12	Mn	−4.62	2.73	2.60	4.05	5.3514

---

(c)
Host	Dp	M	μnet (µB)	Δgxx	Δgyy	Δgzz	geff	μeff (µB)
Ga9S12	Zr	4	3.00	−1.076	−0.998	−0.146	1.26	2.44
Ga9Se12	Mo	6	5.00	−0.434	−0.275	−0.010	1.76	5.22
Ga9Te12	Mo	8	7.00	0.494	1.141	5.862	4.50	17.86
In9S12	Zr	2	1.00	−0.204	−0.204	−0.003	1.87	1.62
In9Se12	Mo	8	7.00	0.048	0.424	0.650	2.38	9.43
In9Te12	Mo	8	7.00	0.884	0.892	8.042	5.28	20.93
Tl9S12	Mn	7	6.00	0.032	0.149	0.231	2.14	7.41
Tl9Se12	Mn	7	6.00	0.947	1.632	2.481	3.70	12.78
Tl9Te12	Mn	5	4.00	−1.232	−1.207	−0.073	1.17	2.85

---

After introducing the TM atoms at the central site, the formation energy is computed as

E_form = E(M₉TMX₁₂) + E_M − E(M₁₀X₁₂) − E_TM, (1)

where E(M₉TMX₁₂) is the total energy of the defected systems, E_M is the total energy of an isolated host metal atom, E(M₁₀X₁₂) is the total energy of the pristine system, and E_TM is the total energy of an isolated TM atom. As depicted in Table 2(b), the derived substitution energies are found to be exothermic in nature and E_form ranges between −4.2 eV and −5.17 eV, i.e., the TM introduction is energetically well favorable. Such energetic favorability is also observed for their corresponding charged states. The corresponding formation energies are provided in Table S1. In addition to deriving the binding and formation energies, the physical stability of the pristine and perturbed configurations was further confirmed via frequency calculations and ab initio molecular dynamics. The vibrational analysis provides IR spectra, where spectral analysis shows that all vibrational frequencies are real (positive), with no imaginary (negative) modes, which indicates that the configurations are dynamically stable and do not exhibit time-divergent atomic displacements (Fig. S3, SI). The later calculations are performed to examine whether the observed dynamic stability is also maintained under thermal perturbation (Fig. S4, SI). The time-evolution of potential energy of all the structures is evaluated over 5 femtoseconds at an instantaneous perturbation temperature of 300 K. The results show a linear relationship between potential energy, with no abrupt large increase or decrease, which confirms the thermal stability of the nanodots at room temperature. These findings demonstrate that the nanodots remain both dynamically and thermally stable over the investigated temperature range. As a result, experimental synthesis of the pristine and substituted configurations appears feasible under conventional laboratory conditions.

As observed from the PDOS analysis in Fig. 3, the low-lying s–p states are shifted and hybridized with the TM d states, and these overlaps are visible around the HOMO and LUMO regions. The orbital overlap increases compared to the pristine's orbital hybridization, and the increased p-orbital contribution mainly originates from the high-energy p orbitals introduced by the Ga/In/Tl atoms. The energy levels are more densely spaced around the energy gap region, particularly for β-channels, which causes the energy gap to be lower in that channel for most of the configurations. The pristine nanodots adopt closed-shell singlet ground states in SCF calculations, with all electrons paired in occupied orbitals. The derived HOMO–LUMO gap ranges from 0.89 to 5.02 eV, according to Table 2. For Ga-based nanodots, the highest energy gap is shown by Ga₁₀Te₁₂ (3.5 eV) and the gap keeps reducing when the chalcogen changes from Te to S. For the In-based QDs, the gap increases when S → Se (by 0.9 eV) but slightly decreases when Se → Te (by 0.2 eV), while for the Tl-based QDs it sharply decreases by 2.17 eV (for S → Se) and then again decreases by 1.55 eV (for Se → Te). Such non-monotonic behavior observed in In- and Tl-based QDs is primarily due to their larger orbital energy (5p¹, 6p¹) compared to the smaller orbital energy of Ga-based MXs (4p¹). The larger orbital tends to be more diffusive, and these diffuse orbitals extend further into space, which can increase the orbital interaction range. The replacement with TM atoms causes the energy gaps to split into spin-up (E^↑_g) and spin-down (E^↓_g) orbital channels, and almost for all the systems the E^↓_g channel is more energetically accessible, as presented in Table 2(b). Also a considerable widening of 4.41 (E^↑_g) and 4.01 (E^↓_g) eV is observed for Ga₉S₁₂ QDs upon Zr substitution. However, this shift notably reduces at E^↑_g (by 1.51 eV) and E^↓_g (by 1.46 eV) when the host metal is replaced by indium, which indicates host-dependent orbital hybridization for the chosen dopants. The dipole moment analysis from Table 1(a) supports such host-dependent orbital hybridization by quantifying the charge redistribution upon incorporation of the dopant. Table 2 shows that inclusion of Zr heavily alters the charge distributions (by nearly three orders of magnitude) in the Ga–S framework in contrast to the In–S framework (only ∼0.67 times). Moreover, a reasonable energy gap increment of 3.06 (E^↑_g), 0.8 (E^↓_g), and 1.83 (E^↑_g), −1.88 (E^↓_g) eV is observed for Ga₉Se₁₂ and Ga₉Te₁₂, respectively, upon introducing a 4d⁵ perturbation. In contrast, the corresponding increments for In₉Se₁₂ and In₉Te₁₂ are relatively modest, at approximately 0.43 (E^↑_g), −3.1 (E^↓_g) eV, and 0.9 (E^↑_g), −2.95 (E^↓_g) eV. A decrease in the E^↑_g-channel is found for both Tl₉S₁₂ and Tl₉Te₁₂ by 0.17 and 1.8 eV, while the E^↓_g-channel drops by 3.08 and 0.35 eV. For Tl₉Se₁₂, E^↑_g increases by 1.81 eV, whereas E^↓_g decreases by 0.97 eV. Such variations show how local TM defect sites can differently affect the spin-resolved electronic structures across various host compositions.

3.2 Frontier orbital characteristics

Here, we attempt to analyze the spin-dependent orbital distribution in real-space, which is observed to be site-specific according to the isosurface plots. We first briefly compare different spatial configurations, and then perform a mathematical correlation analysis to determine whether such novel distinctions are only visually apparent or physically substantive in origin.

As presented in the isosurface plots of HOMO and LUMO for the pristine (Fig. 4) and substituted (Fig. 5) cases, for the pristine metal-chalcogenide dots, the occupied states are localized at the edge sites as a contribution of the more effective local bonding environment, in contrast to the unoccupied states, which tend to be more dispersed and pronounced around the center region. This particular arrangement results from the extended nature of conduction states that are less tightly bound at these atomic sites. As a consequence, substitution with an unpaired d-orbital atom at this site increases the overall binding energy for the d-orbital-induced hybridization, as depicted in Table 2. According to PDOS analysis, it is also observed that due to such site-specific defect, few spin-dependent localized states appeared within and near the HOMO–LUMO region that has very short amplitudes. These states hybridize with the TM's d-orbitals and form localized wavefunctions that are visible as spatially distinct regions in their corresponding isosurfaces. As a result, these localized regions exhibit not only spin-selectivity but also site-selectivity in physical space. Such spatial-duality is predominantly noticeable for [Ga/In]₉ZrS₁₂, [Ga/In]₉MoSe₁₂, α-Ga₉MoTe₁₂, Tl₉MnS₁₂, α-Tl₉MnSe₁₂, and Tl₉MnTe₁₂ nanoflakes (Fig. 5(a)–(d), (e) and (g)–(i)). Such presence of asymmetric hybridization also arises from the large local dipolar differences induced by TM substitution, as reported in Table 2. We find that the charge redistributions in the aforementioned nanodots are significantly increased relative to their pristine counterparts, with reported values of 1633, 1380, 5254, and 1888 (a.u.) times for selected systems. These outcomes suggest that the coexistence of such varying spin gaps and spatially non-uniform orbitals may allow the nanodots to host multiple energetically distinct electronic states.

To understand whether these molecular dots can represent multiple independent electronic states simultaneously, one way is to verify this theoretically by measuring the spatial correlation between same-spin probability densities. The normalized density–density overlap measures whether spin-dependent electron densities share the same spatial region or are localized at different sites. We can compute this pair-density overlap as

which is normalized by in order to remove the dependence on the absolute spatial localization of individual orbitals and to obtain a bounded, dimensionless measure of relative spatial correlation. As presented in Table 3, the normalized density correlation (D̃_HL) tends to approach zero for most of the QDs. The measured D̃_HL confirms our intuition that two electron-density distributions that are separated by a finite energy gap are apparently orthogonal (essentially ≤10⁻³) in terms of their real-space envelopes for most of the dots. However, this correlation measurement of identical spin-densities with |ψ|² can only capture the positive phase and disregards signed-amplitudes or nodal structures.

Table 3 Spatial density–density correlation between same-spin orbitals of different energies

System	Channel	Overlap pair	DHL	D̃HL
Ga9ZrS12	α	Hα–Lα	3.3888 × 10^−4	0.0746
	β	Hβ–Lβ	3.1734 × 10^−5	0.0103
Ga9MoSe12	α	Hα–Lα	5.1302 × 10^−5	0.0104
	β	Hβ–Lβ	3.2719 × 10^−5	0.0123
Ga9MoTe12	α	Hα–Lα	1.5198 × 10^−4	0.0594
	β	Hβ–Lβ	8.9530 × 10^−4	0.7851
In9ZrS12	α	Hα–Lα	1.0305 × 10^−4	0.0350
	β	Hβ–Lβ	1.5528 × 10^−5	0.0043
In9MoSe12	α	Hα–Lα	2.4522 × 10^−4	0.0536
	β	Hβ–Lβ	1.4173 × 10^−3	0.9492
In9MoTe12	α	Hα–Lα	2.2926 × 10^−4	0.1213
	β	Hβ–Lβ	7.4158 × 10^−4	0.7850
Tl9MnS12	α	Hα–Lα	1.0665 × 10^−4	0.0489
	β	Hβ–Lβ	3.1691 × 10^−5	0.0109
Tl9MnSe12	α	Hα–Lα	1.6890 × 10^−4	0.0928
	β	Hβ–Lβ	1.7834 × 10^−3	0.8797
Tl9MnTe12	α	Hα–Lα	1.8716 × 10^−5	0.0081
	β	Hβ–Lβ	6.1843 × 10^−5	0.0477

---

This leads us to measure another quantity which can characterize the real-space interrelation of spatial phases and nodal structures between two same-spin orbitals, and is rotation-invariant. Because two orbitals may occupy the same spatial region (as per the overlapping envelopes) but have opposite signs or nodes locally, it may be possible to address them selectively or separately (as they correspond to energetically distinct levels). In order to verify this, the spatial overlap between same-spin orbital wavefunctions is evaluated from the reconstructed three-dimensional wavefunction grids ψ(x,y,z), that are derived from the signed real-space orbital amplitudes ψ(r). The overlap integral S quantifies the spatial correlation between HOMO and LUMO orbitals within a spin channel (which is basically 〈ψ_H|ψ_L〉), where the normalized overlap S_norm accounts for orbital scaling. Within each spin channel, the overlap integral

is computed on the common cubic grid, and normalization by
yields S_norm. A value of |S_norm| < 10⁻⁴ (that is near zero) would indicate that the two same-spin orbitals are spatially orthogonal. A normalized overlap approaching near zero indicates that the same-spin orbitals that are separated by a well-defined energy gap do not result in significant spatial mixing or hybridization in space (i.e., the asymmetry is distinguishable). As depicted in Table 4, for the same-spin orbital pairs the normalized overlap |S_norm| < 10⁻⁴ which confirms their spatial orthogonality. This behavior also aligns with the presented isosurface plots, which visually confirm the non-overlapping and site-selective nature of these orbital shapes depending only on their spin.

Table 4 HOMO–LUMO overlap integrals from signed real-space wavefunction, normalized overlaps, and electronic spatial extents (〈r²〉 in Bohr²) for all structures and spin channels

System	Channel	Overlap pair	S	|Snorm|	〈r^2〉H	〈r^2〉L
Ga9ZrS12	α	Hα–Lα	1.2474 × 10^−6	0.0000	26.5133	73.4516
	β	Hβ–Lβ	1.7430 × 10^−6	0.0000	182.1451	165.0047
Ga9MoSe12	α	Hα–Lα	−3.8323 × 10^−5	0.0000	188.4636	30.0014
	β	Hβ–Lβ	−1.9331 × 10^−6	0.0000	188.3397	184.7601
Ga9MoTe12	α	Hα–Lα	−1.2736 × 10^−6	0.0000	80.1892	99.1496
	β	Hβ–Lβ	−9.4803 × 10^−6	0.0000	189.8775	202.7604
In9ZrS12	α	Hα–Lα	2.6031 × 10^−5	0.0000	215.4526	57.3589
	β	Hβ–Lβ	−8.6911 × 10^−6	0.0000	215.4633	18.0892
In9MoSe12	α	Hα–Lα	−1.3074 × 10^−4	0.0001	33.8344	94.8511
	β	Hβ–Lβ	2.2814 × 10^−4	0.0002	203.5608	214.1154
In9MoTe12	α	Hα–Lα	−7.8615 × 10^−5	0.0001	99.3915	114.9374
	β	Hβ–Lβ	−6.7803 × 10^−4	0.0007	214.9433	228.5232
Tl9MnS12	α	Hα–Lα	6.4360 × 10^−6	0.0000	192.6453	99.3229
	β	Hβ–Lβ	1.1188 × 10^−6	0.0000	193.1523	191.1198
Tl9MnSe12	α	Hα–Lα	1.7724 × 10^−4	0.0002	206.9119	111.9279
	β	Hβ–Lβ	−1.4000 × 10^−5	0.0000	200.9388	198.8377
Tl9MnTe12	α	Hα–Lα	3.4190 × 10^−5	0.0000	228.2090	182.6864
	β	Hβ–Lβ	6.4570 × 10^−7	0.0000	228.2522	146.3248

---

Such existence of spin-selective orbital reconfiguration can further be understood if we observe how the atomic spins are spatially populated, as depicted with the spin-density plot (Fig. 6). It is observed that due to the presence of nested TM–X₃ hybridization in the host configuration, the distribution of spin-up (α) and spin-down (β) electronic states varies significantly. This implies that the physical occupation of different spin-dependent orbitals in real-space would be vastly different, and as we can see from the isosurface plots it is entirely selective of the host structure. Such variations of spin-distribution provide a core intuition that if we introduce site-specific defects with a particular dopant (that has unpaired d-electrons) in the nonmagnetic host to form localized hybridization, it creates a major imbalance of spins; hence it is possible to achieve mapping of electronic orbitals in a way that they would distribute asymmetrically and disparately in real-space. This intuition is consistent with our findings, as not only spin-density plots but also the PDOS analysis suggests (discussed earlier) that the number of different orbital and atomic electronic states are distributed at distinct energy levels due to the presence of Zr–S₃, Mo–(Se/Te)₃, and Mn–(S/Se/Te)₃ hybridization. As a consequence, the α-spin and β-spin orbitals occupy different spatial regions while still hybridizing energetically, and exhibit site-selective localization, for which the reconfigurations are real, rational, and spin-dependent. The spin-selective orbital reconfiguration observed in these nanodots originates from a combined effect of the dopant electronic structure and local hybridization geometry. First, the d-orbital occupancy of the TMs (3d⁵, 4d², 4d⁵) determines how many spin-polarized states are available to participate in bonding, which governs whether the resulting orbital reconfiguration is strong, weak, or absent. This electronic configuration gives rise to an exchange splitting between spin-up and spin-down components, and the magnitude of such splitting then energetically separates the E^↑_g and E^↓_g channels and enforces spin selectivity at the level of the frontier orbitals. Moreover, the relative overlap between the p- and d-orbital energy of TM and chalcogens controls which orbitals hybridize and within which energy window, which in turn causes different TM–X₃ motifs to exhibit distinct reconfiguration patterns. It is further evident that the local TM–X₃ coordination symmetry lifts the degeneracy among spatial orbitals, and as a result, the hybridized states are shaped into distinct spatial orientations and nodal structures. Together, these factors give rise to a logical, site-selective, and spin-dependent orbital reconfiguration, which is consistently supported by spin-density distributions, PDOS features, and real-space isosurface plots.

To further support the real-space orbital and spin-density analyses, we performed Mulliken population analysis for all pristine and centrally perturbed nanodots (Tables S4–S12) along with their representative charged state (Tables S13–S21). It is observed that the resulting charge and spin populations are consistent with the site-selective and spin-dependent electronic reconfigurations identified from the orbital isosurfaces and overlap analyses. Because Mulliken populations project the spin-resolved density matrix on to atom-centered basis functions, the observed atom-resolved spin moments represent the total integrated consequence of the spin-selective orbital reconfigurations, which are visually observed in real-space. In the subsequent analysis, we discuss that, apart from the orthogonality of spin-dependent energy states, the same-spin orbitals also encompass distinct spatial radius, which would provide a more constructive understanding of their spatially independent nature. Collectively, these results imply that most of the nanodots demonstrate quasi-independent like electronic states within their same-spin channel, which could inform future experimental efforts.

3.3 Orbital spatial extent

The calculation of the electronic spatial extent (ESE) provides an outlook into the spatial characteristics of orthogonal wavefunctions, enabling the quantification of how each orbital's electron density is distributed in space. It is the expectation value of the squared radial distance 〈r²〉, which is derived from the orbital density aswhere ΔV is the voxel volume and the summation over all voxels approximates the integral over space. The continuous integral is then discretized on the 3D cube grid by multiplying each voxel's contribution by its volume and summing over all voxels, which were normalized by the total orbital density to account for non-normalized wavefunctions.

In the case of different energy orders with the same spins, the quadrupole moment analysis can give a more constructive outlook (see Table S2 in the SI), as it quantifies the degree of spatial asymmetry and directional electron density distortion across the molecular framework, and indicates a more spatial polarization and anisotropy in the substituted configurations. As per Table 4, for the same energy order with different spins, we can see that due to the Zr–S₃ interaction within the host of Ga₉S₁₂ QDs, ESE shifts about +155.6 and +91.5 a.u. for the spin-up and spin-down orbitals, respectively. However for the same hybridization, In₉S₁₂ shows an identical distribution of E^↑_g electrons and a reduction of 39 a.u. for E^↓_g. Moreover, configurations like α-Ga₉MoS₁₂, α/β-In₉ZrS₁₂, α-Tl₉MnS₁₂, and α-Tl₉MnSe₁₂ have a comparatively large spatial difference of 158.4, 158.1, 197.37, 93.3, and 95 (a.u.), respectively, which tells us that alongside their orthogonal behavior the expansion in real-space deviates such that each of the orbitals can possibly be realized as a quasi-independent state (in real-space).

For the rest of the dots, the ESE shifts positively and Table S3 (in the SI) verifies whether the spatial distribution of both spins within the same energy window is identical. We find that the occupied orbitals of Ga₉MoSe₁₂, In₉ZrS₁₂, Tl₉MnS₁₂, and Tl₉MnTe₁₂ are non-orthogonal in nature and normalized as

Despite such similar distributions, the ESE changes significantly for the virtual orbitals of both spins. This reduction or increment can be understandable if we again observe the spatial dispersion of spin-differential frontier orbitals (as per Fig. 5). For example, in the case of the Ga₉MoSe₁₂ nanodots, we see that HOMOs are distributed along the same spatial location, simultaneously for both channels, but as the dipole moments are increased by 1630 times from the pristine, the spatial extent of electrons increases by some amount (+155 a.u. at LUMO in this case). But this explanation may seem contradictory at first look if we observe In₉ZrS₁₂, because ESE shifted negatively despite its dipole moment remained the same compared to the pristine. The main reason this result is not counterintuitive is that the local electronic stabilization (dopant-induced potential wells that preferentially localize electrons in certain regions) has more effect on orbital localization than the global dipole variation. Moreover, though Ga₉MoSe₁₂, In₉ZrS₁₂, Tl₉MnS₁₂, and Tl₉MnTe₁₂ exhibit a clear homogeneity in the spatial localization of α-HOMO states localized below the orbital energy gap upon substitution (Fig. 5(b), (d), (g) and (i)), they confirmed higher spin-states like quintet, sextet, and septet. Notably, distinct spin-dependent energy gaps are still present despite this spatial homogeneity. The energy gap difference between alpha and beta channels is 2.26, 0.05, 2.91, and 1.45 (eV), respectively, as reported in Table 2.

Such coexistence of a small but finite density overlap with a vanishing signed wavefunction overlap, as well as distinct spatial extent, shows that these states remain functionally decoupled due to site- and node-selective separation, rather than only phase cancellation within the same region. Moreover, though PDOS depicts the presence of energetic hybridization between dopant and host orbitals, the overlap results show that these states remain spatially distinct and occupy different sites and spatial radii. This spatial separation limits real-space coexistence despite energy-space mixing, which in turn allows the spin-dependent orbitals to behave as distinguishable spatial channels. Such real-space localization suggests that selective excitation or manipulation of site- and spin-specific states could be feasible, which holds promise for memory or logic functionalities at such nanoscale limits. Taken together, the observed anisotropic spatial distribution of orbitals allows us to control each of the parameters independently and suggests that such spin-dependent, site-specific confinement provokes differential localization of orbital polarization, which in turn leads to selective spatial distribution of frontier orbitals and spin channels.

3.4 Magnetic strength and spin-dependent response

In this subsection, the goal of the discussion is to structure a constructive analysis on the observed magnetic response, and how the intrinsic as well as the field-induced magnetic behavior recharacterizes spin-dependent interactions in the system. The spin multiplicity (M) is obtained by performing the ground state geometry optimization for different spin-polarized states (quartet, septet, octet, etc.), and the one with the lowest total energy is defined as the magnetic state of that structure. For all configurations studied, we explicitly tested multiple spin multiplicities. We considered states with up to 16 unpaired electrons. In every case, the total energy shows that the highest-spin state is the most stable. This ensures that the identified ground state indeed corresponds to the highest-spin configuration. In all cases, the computed total energies confirmed that the highest-spin state is the most stable, which validates our assumption. The net magnetic moments μ_net are derived from the Mulliken spin analysis which represents the sum of all Mulliken spins induced in the system. As depicted, the TM-dopant carries the maximum-site-moment and its magnitude depends on the host elements. As we can see from Tables S4–S12 (SI), the Zr–S₃ interactions in Ga₉S₁₂ produce a maximum local moment of 1.935µ_B, and when the host framework constitutes In₉S₁₂ it shows 1.048µ_B. But for the Mn–S₃ interactions in Tl₉S₁₂ the source magnetic moment can reach up to 4.742µ_B. This is because when the orbitals are more spatially extended, the spatial overlap with dopant d-orbitals becomes more frequent and higher, which in turn leads to more hybridization per site. Such observation is similar for the other nested hybridization, as we can see that Mn–Se₃ and Mn–Te₃ tend to produce comparatively higher or equivalent local moments compared to the Mo–Se₃ and Mo–Te₃ interactions with the host.

The magnetic responses of the nanodots when interacting with an external magnetic field are presented in Table 2(c). As presented, the significant shifts in directional g-tensors which increase the effective g-factor are the direct consequences of the high spin–orbit coupling (SOC) effect in the nanoflakes resulting from the combined orbital's contribution from transition metals and heavy metals to magnetism being much higher and dominating over the spin's contributions. In response to the interaction with a magnetic field, the QDs exhibit an effective magnetic moment which can be quantified as . g_eff is derived from the g-shifts of each principal axis and calculated as their averaged value, where μ_eff increases along with the increase in g_eff. The computed g-tensors include contributions from orbital-Zeeman and perturbative spin–orbit terms (g_OZ/SOC), and the observed anisotropic g-factor shifts primarily arise from orbital-Zeeman contributions, where SOC provides a quantitative correction. Deviations from the free-electron g-factor indicate enhanced orbital participation and anisotropic magnetic response induced by the dopants. The mean g-factor of the free electron is a key descriptor to understand the spin behavior, where g_e ∼ 2.002 and any deviation from this value indicates modifications of SOC and exchange interaction in the molecular electronic configurations due to central defects. The deviation of g-factor along the principal axis is derived with the GIAO formalism, and the details are presented in the methodology section of the SI. As we can see, when the effective g-factor increases, μ_eff exhibits higher response, which reflects its ability to enhance orbital-induced magnetism due to the strong presence of SOC. The field-induced magnetic moments due to such effect show very strong magnetic moments, as the derived μ_eff ranges between 1.62µ_B and 20.93µ_B. This kind of large magnetic response suggests multifunctional magnetic behaviors of these nanoflakes, which show the simultaneous presence of spin-selective orbital distribution in real-space and field-controllable magnetic moments.

The PDOS for all the considered quantum dots is shown in Fig. 3a–c, which depict the specific orbital contributions involved in host-TM hybridization. For the pristine gallium- and indium-based MX dots, it is observed that metal s, p and chalcogen p orbitals dominate primarily because the d-orbitals remain stable (appear at very low energies) and do not participate in bonding, staying inactive near or around the energy gap region. Compared to the thallium-based QDs, electrons from the d-orbitals actively participate in the PDOS spectra of the pristine system. This behavior arises from relativistic effects, particularly scalar relativistic contraction of the 6s orbital and expansion of the 5d shell, which raise the energy of the 5d states and enhance their spatial extension. As a result, Tl's 5d orbitals can hybridize appreciably with neighboring p orbitals (from chalcogen atoms and the host configuration), and show a detectable d-character near the HOMO–LUMO region of the electronic spectrum. A tiny contribution comes up from Zr's 4d² orbitals in the Ga₉S₁₂ and In₉S₁₂ nanodots upon single Zr substitution (Fig. 3a), showing a marginal overlap with the host orbitals. It is worth noting that this limited interaction still contributes to the spin-polarization observed in the quartet (M = 4) and doublet (M = 2) magnetic states (spin-density isosurface of Fig. 6a and d), particularly for the In₉S₁₂ system, where the spin density is closely confined at the substituted site. In contrast to the Tl₉S₁₂ QDs, Mn substitution results in more noticeable orbital overlap due to its higher number of unpaired d-electrons (3d⁵), and in both spin channels this interaction is visible along the HOMO–LUMO (Fig. 3c) region. Despite the lower availability of electronic states, the interactions are effective enough to produce high septet polarization, and the resultant magnetic moments show a broad spatial distribution from the center to an in-plane edge site adjacent to the substituted site (Fig. 6g).

Following Mo-induced perturbation in Ga₉Se₁₂, Ga₉Te₁₂, In₉Se₁₂, and In₉Te₁₂ nanoflakes, the orbital interactions lead to much higher spin-polarized states. The presence of 4d⁵ orbitals, which hybridize in the midgap region of Ga₉Se₁₂ (Fig. 3a), results in asymmetric spin localization around the dopant and at an edge site (Fig. 6b), similar to what is observed for Ga₉S₁₂ and Tl₉S₁₂ nanoflakes. This sort of spatial separation of spin densities at a peripheral site on the same plane as the dopant is evident in Tl₉Se₁₂ (Fig. 6h) and Tl₉Te₁₂ (Fig. 6i) QDs, where d-orbital-induced hybridization via Mn substitution is noticeably stronger and denser (Fig. 3c).

But a more widespread delocalization of magnetic moment, which extends from the central to peripheral edge sites across the plane, is observed for Ga₉Te₁₂ (Fig. 6c), In₉Se₁₂ (Fig. 6e), and In₉Te₁₂ (Fig. 6f) nanoflakes. The nested Mo–(Se, Te)₃ interactions shift the low-lying valence states upward in energy through 4d⁵ orbital overlap (Fig. 3a and b), and consequently, the magnetic configuration corresponds to the octet state, which is the highest among the group (Table 2). It is understandable from our PDOS analysis that when the number of overlapping states from hybridized orbitals increases near or above the energy gap, exchange interactions become more strong, which in turn leads to the emergence of site-specific spin textures and higher spin polarization in the system. To provide a more detailed outlook on atomic contributions, the element-projected PDOS is presented in the SI (Fig. S1).

3.5 Static optical response

The anisotropy in α(0;0,0), β(0;0,0), and γ(0;0,0) has been studied in response to the TM defect in order to characterize the nature of the material and to understand whether the observed anisotropic magneto-electronic response is also present when interacting with a constant electric field, at such confined geometry. The ground-state electron cloud response increases when the chalcogen changes from S to Se, Te, as measured by the average linear polarizability (α_iso) and directional polarizability (α_aniso) for pristine QDs. The relative ratio between anisotropic and isotropic responses (R_A) provides a quantitative measure of the directional dependency in these nanoflakes. The isotropic and anisotropic responses range from ∼0.6 × 10³ to ∼0.15 × 10⁴ for the pristine configurations, with values on the order of ∼10⁴ observed when the MX dots are bonded with Te, as depicted in Table 5. This behavior is a direct consequence of shifted dipole moments that increase when either the metal or chalcogen atom is replaced by a higher atomic radius atom. Compared to the substituted configurations, the anisotropic response shows a marginal increase for Ga₉MoTe₁₂, In₉MoSe₁₂, In₉MoTe₁₂, Tl₉MnS₁₂, and Tl₉MnTe₁₂, but significantly dominates for Tl₉MnSe₁₂, reaching up to 0.22 × 10⁴ scale.

Table 5 Calculated static linear polarizability and higher-order hyperpolarizabilities (α, β, γ) for the considered quantum dots

Host	Dp	αiso (a.u.)	αaniso (a.u.)	β‖ (a.u.)	β⊥ (a.u.)	γ‖ (a.u.)	γ⊥ (a.u.)
Ga10S12	—	6.1 × 10^2	6.0 × 10^2	−3.8 × 10^1	−1.3 × 10^1	3.7 × 10^5	1.2 × 10^5
Ga9S12	Zr	6.1 × 10^2	5.6 × 10^2	1.8 × 10^3	6.1 × 10^2	6.1 × 10^5	2.0 × 10^5
Ga10Se12	—	7.4 × 10^2	7.7 × 10^2	5.6	1.9	4.3 × 10^5	1.5 × 10^5
Ga9Se12	Mo	7.5 × 10^2	5.0 × 10^2	4.5 × 10^3	1.5 × 10^3	8.2 × 10^5	2.7 × 10^5
Ga10Te12	—	1.0 × 10^3	1.1 × 10^3	2.8	9.4 × 10^−1	6.2 × 10^5	2.1 × 10^5
Ga9Te12	Mo	1.1 × 10^3	1.3 × 10^3	2.2 × 10^5	7.2 × 10^4	3.7 × 10^8	1.2 × 10^8
In10S12	—	7.3 × 10^2	7.1 × 10^2	3.5 × 10^1	1.2 × 10^1	9.1 × 10^5	3.1 × 10^5
In9S12	Zr	7.5 × 10^2	6.9 × 10^2	2.0 × 10^1	6.5	9.3 × 10^5	3.1 × 10^5
In10Se12	—	8.8 × 10^2	8.8 × 10^2	−1.5 × 10^1	−5.0	2.7 × 10^5	8.0 × 10^5
In9Se12	Mo	9.2 × 10^2	9.7 × 10^2	1.4 × 10^5	4.7 × 10^4	9.9 × 10^7	3.3 × 10^7
In10Te12	—	1.1 × 10^3	1.2 × 10^3	−7.6 × 10^1	−2.5 × 10^1	1.1 × 10^6	3.5 × 10^5
In9Te12	Mo	1.2 × 10^3	1.3 × 10^3	2.2 × 10^5	7.3 × 10^4	1.6 × 10^8	5.2 × 10^7
Tl10S12	—	8.4 × 10^2	9.0 × 10^2	2.7 × 10^1	9.2	8.0 × 10^5	2.7 × 10^5
Tl9S12	Mn	8.2 × 10^2	8.7 × 10^2	3.0 × 10^3	1.0 × 10^3	9.9 × 10^5	3.3 × 10^5
Tl10Se12	—	1.0 × 10^3	1.1 × 10^3	2.2 × 10^1	7.2	8.2 × 10^5	2.7 × 10^5
Tl9Se12	Mn	1.7 × 10^2	2.2 × 10^3	−5.6 × 10^4	−1.9 × 10^4	−1.2 × 10^8	−3.9 × 10^7
Tl10Te12	—	1.3 × 10^3	1.5 × 10^3	1.2 × 10^2	4.0 × 10^1	1.1 × 10^6	3.8 × 10^5
Tl9Te12	Mn	1.3 × 10^3	1.4 × 10^3	2.3 × 10^3	7.6 × 10^2	2.3 × 10^5	7.6 × 10^4

---

This large increment due to Mn–Se₃ can be rationalized upon inspection of the corresponding isosurface plots of the HOMO (Fig. 5(c), (e), (f), (g) and (i)). They show a highly asymmetric localization of wavefunctions and diffuse orbitals, where the spin-up orbitals are primarily confined to one edge site, while the spin-down orbitals are distributed along a different peripheral site in the opposite direction (yet along the same axis). Such directional separation from the overlap of 6p¹–3d⁵–4p⁴ intensifies charge polarization along one axis, which in turn amplifies the anisotropic response of this quantum dot.

More interestingly, the first-hyperpolarizability tensor components (β_‖, β_⊥) exhibit phase (direction of the induced polarization) inversion between the applied electric field and the induced polarization. As shown in Table 4, this phase inversion becomes switchable to in-phase polarization through Zr and Mo substitution in pristine Ga₁₀S₁₂, In₁₀Se₁₂, and In₁₀Te₁₂, while Mn substitution in Tl₁₀Se₁₂ induces out-of-phase polarization relative to its pristine counterpart. The in-phase polarization increased by roughly ∼10² to ∼10⁴ times of magnitude, while the out-of-phase response reached up to ∼10³ times of magnitude in their corresponding β_‖ components. Alongside this, the Ga₉MoTe₁₂ nanodot exhibits the highest increase, increasing approximately six orders of magnitude, and the largest variation between the parallel and perpendicular components is observed for In₉ZrS₁₂ (R_B reaches about 3.06).

In addition, a similar phase reversion in second-hyperpolarizability tensor components (γ_‖, γ_⊥) is also observed for the Tl₉MnSe₁₂ nanodot, where the intensity reaches approximately 10⁹ and 10⁸ along the field direction and perpendicular to it, respectively. As reported in Table 5, the third-order nonlinear response is the strongest among all considered QDs, which suggests that such confined geometry, zigzag edge-bonding, and local perturbations give rise to long-range, non-saturating, and significant static nonlinearity in these sub-de Broglie length regimes. The highest anisotropy ratio (R_Γ) is observed for the Ga₁₀S₁₂ nanoflake (∼3.1), which indicates that its wavefunction delocalization (Fig. 4(a) and (b)) is much more pronounced along the direction of the interacting field. The directional uniformity preserved in the systems, in accordance with such amplified intensities, can be more easily understood from Fig. 7.

The predicted spin- and site-selective orbital reconfigurations could be experimentally validated in several ways.^49–52 Spin-polarized scanning tunneling microscopy (SP-STM) can directly visualize the real-space spin densities and reveal the spatial separation between E^↑_g and E^↓_g orbitals. Another method could be X-ray magnetic circular dichroism (XMCD), which can probe the spin-dependent contributions of the transition metal dopants and confirm magnetic anisotropy. Moreover, spin-resolved photoemission spectroscopy (spin-PES) can measure the energy-level separation between spin-up and spin-down orbitals, which may directly validate the predicted exchange splitting. In addition to this, STM or scanning tunneling spectroscopy (STS) can map the local density of states to confirm site-specific localization of frontier orbitals. X-ray absorption spectroscopy (XAS) or electron energy-loss spectroscopy (EELS) can also probe the symmetry and spatial distribution of unoccupied orbitals. The effective g-factor enhancement due to the SOC effect can be confirmed via electron spin resonance (ESR) experiments under applied magnetic fields. Furthermore, linear polarizability and first/second-hyperpolarizabilities can be measured using techniques like electric-field-induced second-harmonic generation (EFISHG), optical absorption spectroscopy, or hyper-Rayleigh scattering (HRS) measurements, which provide information on charge distribution and anisotropy in orbital delocalization. Taken together, these findings highlight that engineered point-defect perturbations lead to systematic tuning of the magnitude, directionality, and phase of static linear and nonlinear optical responses through controlled spin–orbital reconfiguration, where spin-selective orbital localization, edge-confined geometry, and TM-dopant identity collectively govern these properties across sub-de Broglie length scale quantum dots.

4 Conclusion

Using DFT, the spatial regulation of spin-dependent frontier orbitals and anisotropy in hyper(polarizabilities) of group-13 based one-atom-thick zigzag-edged PTMC quantum dots is investigated. A metal atom from the host system was substituted with a single TM atom, which forms a nested TM–X₃ region, and a systematic computational approach was used to determine how such confined interactions within a closed-shell structure modify the morphology of electronic orbitals and thus their spatial distributions, spin decoration, and optoelectronic response. Spin-dependent spatial expansion and contraction of wavefunctions were observed, along with edge-centric and spin-down dominated magnetic density and the dominance of longitudinal or transverse hyperpolarizabilities. This study establishes that point-defect perturbations in atomically thin dots can introduce quasi-independent electronic states in real-space and induce phase inversion of nonlinear polarization, which enables control over diverse multifunctional properties important for multivariable quantum information processing technologies.

The edge atoms of the QDs were modeled in their under-coordinated form (without applying explicit passivation) to isolate the intrinsic effects of internal substitutional defects, while TM dopants were placed at a central interior site to ensure a controlled coordination environment. Future work could explore the effects of alternative defect configurations or host stoichiometry, which may include different TM incorporation sites, vacancies, and edge substitutions with appropriate chemical passivation, to establish broader structure–property trends in QDs. Such studies would extend insights into how defect configurations influence the electronic and optical behavior of quantum dots.

Conflicts of interest

The authors declare that they have no known competing financial or personal interests that could have appeared to influence the work reported in this paper.

Data availability

Supplementary information: additional details and data to complement the main text, including computational methods, structural parameters, formation energies of charged states, charge distribution analyses, spin and orbital projections, total density of states, atomic charges and spin populations, vibrational responses, finite-temperature energy landscapes, and full atomic coordinates of the studied quantum dots. See DOI: https://doi.org/10.1039/d6ma00156d.

Information on supporting data from this study remains privately accessible to the authors, who can provide it to qualified requesters upon request.

Acknowledgements

The authors declare that this research was conducted without any specific financial support from public, commercial, or not-for-profit funding agencies.

References

1. B. Dieny, I. L. Prejbeanu, K. Garello, P. Gambardella, P. Freitas, R. Lehndorff, W. Raberg, U. Ebels, S. O. Demokritov, J. Akerman, A. Deac, P. Pirro, C. Adelmann, A. Anane, A. V. Chumak, A. Hirohata, S. Mangin, S. O. Valenzuela, M. C. Onbaşli, M. D’Aquino, G. Prenat, G. Finocchio, L. Lopez-Diaz, R. Chantrell, O. Chubykalo-Fesenko and P. Bortolotti, Opportunities and challenges for spintronics in the microelectronics industry, Nat. Electron., 2020, 3, 446–459,  DOI:10.1038/s41928-020-0461-5.
2. A. M. J. Zwerver, T. Krähenmann, T. F. Watson, L. Lampert, H. C. George, R. Pillarisetty, S. A. Bojarski, P. Amin, S. V. Amitonov, J. M. Boter, R. Caudillo, D. Correas-Serrano, J. P. Dehollain, G. Droulers, E. M. Henry, R. Kotlyar, M. Lodari, F. Lüthi, D. J. Michalak, B. K. Mueller, S. Neyens, J. Roberts, N. Samkharadze, G. Zheng, O. K. Zietz, G. Scappucci, M. K. Veldhorst, L. M. K. Vandersypen and J. S. Clarke, Qubits made by advanced semiconductor manufacturing, Nat. Electron., 2022, 5, 184–190,  DOI:10.1038/s41928-022-00727-9.
3. X. Gu, L. Guo, Y. Qin, T. Yang, K. Meng, S. Hu and X. Sun, Challenges and prospects of molecular spintronics, Precis. Chem., 2023, 2, 1–13,  DOI:10.1021/prechem.3c00071.
4. Z. Zhuo, X. Wu and J. Yang, Two-dimensional phosphorus porous polymorphs with tunable band gaps, J. Am. Chem. Soc., 2016, 138, 7091–7098,  DOI:10.1021/jacs.6b02964.
5. A. Ramasubramaniam, D. Naveh and E. Towe, Tunable band gaps in bilayer graphene-BN heterostructures, Nano Lett., 2011, 11, 1070–1075,  DOI:10.1021/nl1039499.
6. Y. Mi, Z. Zhang, L. Zhao, S. Zhang, J. Chen, Q. Ji, J. Shi, X. Zhou, R. Wang, J. Shi Gespräch, W. Du, Z. Wu, X. Qiu, Q. Zhang, Y. Zhang and X. Liu, Tuning excitonic properties of monolayer MoS2 with microsphere cavity by high-throughput chemical vapor deposition method, Small, 2017, 13, 1701694,  DOI:10.1002/smll.201701694.
7. H. Abdelsalam and Q. F. Zhang, Properties and applications of quantum dots derived from two-dimensional materials, Adv. Phys.: X, 2022, 7, 2048966,  DOI:10.1080/23746149.2022.2048966.
8. S. Golovynskyi, M. Bosi, L. Seravalli and B. Li, MoS2 two-dimensional quantum dots with weak lateral quantum confinement: intense exciton and trion photoluminescence, Surf. Interfaces, 2021, 23, 100909,  DOI:10.1016/j.surfin.2020.100909.
9. K. Whitham, J. Yang, B. H. Savitzky, L. F. Kourkoutis, F. Wise and T. Hanrath, Charge transport and localization in atomically coherent quantum dot solids, Nat. Mater., 2016, 15, 557–563,  DOI:10.1038/nmat4576.
10. A. Pal, S. Zhang, T. Chavan, K. Agashiwala, C. Yeh, W. Cao and K. Banerjee, Quantum-engineered devices based on 2D materials for next-generation information processing and storage, Adv. Mater., 2022, 35, 2109894,  DOI:10.1002/adma.202109894.
11. F. Jing, Z. Zhang, G. Qin, G. Luo, G. Cao, H. Li, X. Song and G. Guo, Gate-controlled quantum dots based on 2D materials, Adv. Quantum Technol., 2022, 5, 2100162,  DOI:10.1002/qute.202100162.
12. X. Liu, X. Qin, X. Li, Z. Ding, X. Li, W. Hu and J. Yang, Designing two-dimensional versatile room-temperature ferromagnets via assembling large-scale magnetic quantum dots, Nano Lett., 2021, 21, 9816–9823,  DOI:10.1021/acs.nanolett.1c03814.
13. S. K. Matta, C. Tang, A. P. O’Mullane, A. Du and S. P. Russo, Density functional theory study of two-dimensional post-transition metal chalcogenides and halides for interfacial charge transport in perovskite solar cells, ACS Appl. Nano Mater., 2022, 5, 14456–14463,  DOI:10.1021/acsanm.2c02812.
14. B. Guo, B. Liang, J. Zheng, S. Ahmed, S. Krishna, A. Ghosh and J. Campbell, Digital alloy-grown InAs/GaAs short-period superlattices with tunable band gaps for short-wavelength infrared photodetection, ACS Photonics, 2024, 11, 1419–1427,  DOI:10.1021/acsphotonics.3c01268.
15. H. Bergeron, D. Lebedev and M. C. Hersam, Polymorphism in post-dichalcogenide two-dimensional materials, Chem. Rev., 2021, 121, 2713–2775,  DOI:10.1021/acs.chemrev.0c00933.
16. A. Giri, G. Park and U. Jeong, Layer-structured anisotropic metal chalcogenides: recent advances in synthesis, modulation, and applications, Chem. Rev., 2023, 123, 3329–3442,  DOI:10.1021/acs.chemrev.2c00455.
17. Xin Li, et al., Recent progress on surface chemistry II: Property and characterization, Chin. Chem. Lett., 2025, 36(1), 110100,  DOI:10.1016/j.cclet.2024.110100.
18. Xin Li, et al., Recent progress on surface chemistry I: Assembly and reaction, Chin. Chem. Lett., 2024, 35(12), 110055,  DOI:10.1016/j.cclet.2024.110055.
19. Y. Wang and J. Yi Ferromagnetism in two-dimensional materials via doping and defect engineering, Elsevier eBooks, 2020, pp. 95–124 DOI:10.1016/b978-0-08-102154-5.00008-4.
20. Z. Lei, C. I. Sathish, X. Geng, X. Guan, Y. Liu, L. Wang, L. Qiao, A. Vinu and J. Yi, Manipulation of ferromagnetism in intrinsic two-dimensional magnetic and nonmagnetic materials, Matter, 2022, 5, 4212–4273,  DOI:10.1016/j.matt.2022.11.017.
21. Z. Fei, B. Huang, P. Malinowski, W. Wang, T. Song, J. Sanchez, W. Yao, D. Xiao, X. Zhu, A. F. May, W. Wu, D. H. Cobden, J.-H. Chu and X. Xu, Two-dimensional itinerant ferromagnetism in atomically thin Fe3GeTe2, Nat. Mater., 2018, 17, 778–782,  DOI:10.1038/s41563-018-0149-7.
22. Y. Deng, Y. Yu, Y. Song, J. Zhang, N. Z. Wang, Z. Sun, Y. Yi, Y. Z. Wu, S. Wu, J. Zhu, J. Wang, X. H. Chen and Y. Zhang, Gate-tunable room-temperature ferromagnetism in two-dimensional Fe3GeTe2, Nature, 2018, 563, 94–99,  DOI:10.1038/s41586-018-0626-9.
23. J.-Y. You, B. Gu, S. Maekawa and G. Su, Microscopic mechanism of high-temperature ferromagnetism in Fe, Mn, and Cr-doped InSb, InAs, and GaSb magnetic semiconductors, Phys. Rev. B, 2020, 102, 094432,  DOI:10.1103/physrevb.102.094432.
24. S. Goel, L. D. Anh, S. Ohya and M. Tanaka, Ferromagnetic resonance and control of magnetic anisotropy by epitaxial strain in the ferromagnetic semiconductor (Ga0.8,Fe0.2), Phys. Rev. B, 2019, 99, 014431,  DOI:10.1103/physrevb.99.014431.
25. K. Takiguchi, L. D. Anh, T. Chiba, T. Koyama, D. Chiba and M. Tanaka, Giant gate-controlled proximity magnetoresistance in semiconductor-based ferromagnetic-non-magnetic bilayers, Nat. Phys., 2019, 15, 1134–1139,  DOI:10.1038/s41567-019-0621-6.
26. B. Schuler, J.-H. Lee, C. Kastl, K. A. Cochrane, C. T. Chen, S. Refaely-Abramson, S. Yuan, E. Van Veen, R. Roldán, N. J. Borys, R. J. Koch, S. Aloni, A. M. Schwartzberg, D. F. Ogletree, J. B. Neaton and A. Weber-Bargioni, How substitutional point defects in two-dimensional WS₂ induce charge localization, spin-orbit splitting, and strain, ACS Nano, 2019, 13, 10520–10534,  DOI:10.1021/acsnano.9b04611.
27. P. Andrich, C. F. De Las Casas, X. Liu, H. L. Bretscher, J. R. Berman, F. J. Heremans, P. F. Nealey and D. D. Awschalom, Long-range spin wave mediated control of defect qubits in nanodiamonds, npj Quantum Inf., 2017, 3, 29,  DOI:10.1038/s41534-017-0029.
28. V. Ivády, I. A. Abrikosov and A. Gali, First principles calculation of spin-related quantities for point defect qubit research, npj Comput. Mater., 2018, 4, 13,  DOI:10.1038/s41524-018-0132-5.
29. Y. Pan, J. Yang, S. C. Erwin, K. Kanisawa and S. Fölsch, Reconfigurable quantum-dot molecules created by atom manipulation, Phys. Rev. Lett., 2015, 115, 076803,  DOI:10.1103/physrevlett.115.076803.
30. G. Zhao, W. Ma, S. Yu, J. Zhang and K. Wu, Orbital mixing between colloidal quantum dots and surface-bound molecules, J. Phys. Chem. Lett., 2022, 13, 11892–11898,  DOI:10.1021/acs.jpclett.2c03538.
31. R. C. C. Leon, C. H. Yang, J. C. C. Hwang, J. C. Lemyre, T. Tanttu, W. Huang, K. W. Chan, K. Y. Tan, F. E. Hudson, K. M. Itoh, A. Morello, A. Laucht, M. Pioro-Ladrière, A. Saraiva and A. S. Dzurak, Coherent spin control of s-, p-, d- and f-electrons in a silicon quantum dot, Nat. Commun., 2020, 11, 1,  DOI:10.1038/s41467-019-14053-w.
32. A. L. Kaledin, C. L. Hill, T. Lian and D. G. Musaev, Modulating electronic coupling at the quantum dot/molecule interface by wavefunction engineering, J. Chem. Phys., 2019, 150, 124112,  DOI:10.1063/1.5083056.
33. M. Fang and E.-H. Yang, Advances in two-dimensional magnetic semiconductors via substitutional doping of transition metal dichalcogenides, Materials, 2023, 16, 3701,  DOI:10.3390/ma16103701.
34. H. Hong, C. Wu, Z. Zhao, Y. Zuo, J. Wang, C. Liu, J. Zhang, F. Wang, J. Feng, H. Shen, J. Yin, Y. Wu, Y. Zhao, K. Liu, P. Gao, S. Meng, S. Wu, Z. Sun, K. Liu and J. Xiong, Giant enhancement of optical nonlinearity in two-dimensional materials by multiphoton-excitation resonance energy transfer from quantum dots, Nat. Photonics, 2021, 15, 510–515,  DOI:10.1038/s41566-021-00801-2.
35. S. M. Chowdhury, I. Hossain and M. R. Chowdhury, Interstitial TM-P pairing in P₃ -coordinated wide-gap quantum dots: spin-selective insulating states and enhanced hyperpolarizability, Nanoscale, 2026, 18, 2811–2829,  10.1039/D5NR02041G.
36. R. Rejali, D. Coffey, J. Gobeil, J. W. González, F. Delgado and A. F. Otte, Complete reversal of the atomic unquenched orbital moment by a single electron, npj Quantum Mater., 2020, 5, 1,  DOI:10.1038/s41535-020-00262-w.
37. B. M. Wieliczka, A. L. Kaledin, W. E. Buhro and R. A. Loomis, Wave function engineering in CdSe/PbS core/shell quantum dots, ACS Nano, 2018, 12, 5539–5550,  DOI:10.1021/acsnano.8b01248.
38. Z. Xie, T. Zhao, X. Yu and J. Wang, Nonlinear optical properties of 2D materials and their applications, Small, 2024, 20, 11621,  DOI:10.1002/smll.202311621.
39. X. Wang, A. Zhang, S. Qiu, L. Zhao, Y. Chang, H. Liu and J. Gao, Significant enhancement of the second harmonic generation in Janus α-In/GaSe3 structures, Phys. Rev. B, 2025, 111, 155412,  DOI:10.1103/physrevb.111.155412.
40. Q. Wang and H.-L. Zhang, Two-Dimensional Materials for Nonlinear Optics: Fundamentals, Preparation Methods, and Applications, John Wiley & Sons, 2023 DOI:10.1002/9783527838288.
41. A. Taghizadeh, K. S. Thygesen and T. G. Pedersen, Two-dimensional materials with giant optical nonlinearities near the theoretical upper limit, ACS Nano, 2021, 15, 7155–7167,  DOI:10.1021/acsnano.1c00344.
42. R. W. Boyd, Nonlinear Optics, Academic Press, London, 4th edn, 2020 DOI:10.1016/c2015-0-05510-1.
43. M. Z. Alam, I. De Leon and R. W. Boyd, Large optical nonlinearity of indium tin oxide in its epsilon-near-zero region, Science, 2016, 352, 795–797,  DOI:10.1126/science.aae0330.
44. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery, J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, T. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman and D. J. Fox, Gaussian 09, Revision A.02, Gaussian, Inc., Wallingford CT, 2016.
45. J.-D. Chai and M. Head-Gordon, Long-range corrected hybrid density functionals with damped atom-atom dispersion corrections, Phys. Chem. Chem. Phys., 2008, 10(44), 6615,  10.1039/b810189b.
46. P. J. Hay and W. R. Wadt, Ab Initio Effective Core Potentials for Molecular Calculations. Potentials for the Transition Metal Atoms Sc to Hg, J. Chem. Phys., 1985, 82(1), 270–283,  DOI:10.1063/1.448799.
47. Z. Qu and C. Liu, Intrinsic magnetism in edge-reconstructed zigzag graphene nanoribbons, in Graphene Nanoelectronics: Metrology, Synthesis, Properties and Applications, ed. H. Raza, Wiley-VCH, 2013, pp. 9–28 DOI:10.1002/9781118691281.ch2.
48. V. Kalamse, S. Gaikwad and A. Chaudhari, Computational study of 5d transition metal mononitrides and monoborides using density functional method, Bull. Mater. Sci., 2010, 33, 233–238,  DOI:10.1007/s12034-010-0036-6.
49. P. Beck, Tuning the Spin–Orbit Coupling and the Spin Dynamics of Atomic Chains on Superconductors, Diss. Staats-und Universitätsbibliothek Hamburg Carl von Ossietzky, 2022. https://ediss.sub.uni-hamburg.de/handle/ediss/9900.
50. S. Majumdar, P. Avramov and S. Sakai, Manipulating spins at molecular level: an insight into the ferromagnet-organic interface, World Sci. Reference Spin Organics, 2018, 2, 1–61.
51. B. Wang, Advanced Magnetic Characterization using Electron Microscopy and its Application on Spintronic Devices, PhD thesis, The Ohio State University, 2022. https://www.proquest.com/openview/00ec58e48dff1789bef2f890a1f55d56/1?pq-origsite=gscholar&cbl=18750&diss=y.
52. R. Streubel, et al., Magnetism in curved geometries, J. Phys. D: Appl. Phys., 2016, 49(36), 363001,  DOI:10.1088/0022-3727/49/36/363001.

---