Phase engineering of relaxor ferroelectricity in van der Waals crystal

Domain: neuroscience
Published: 2026-02-09
Editorial outcome: major_revision

Abstract
Phase engineering of two-dimensional van der Waals materials offer a rich landscape to uncover the potential ferroelectric properties. The precisely controlling phases in these reduced dimensional materials with the unique atomic arrangement are crucial. Herein, through a novel phase engineering approach, we synthesize CuInP2(S1-xSex)6 crystal with controllable coexistence of monoclinic and trigonal phases, and experimentally realize relaxor ferroelectricity in low-dimensional materials. Composition tailoring creates an intermediate phase transition region, which facilitates the evolution of ferroelectricity into relaxor ferroelectricity and eventually superparaelectricity. The relaxor ferroelectric properties are verified by the hysteresis-free polarization switching and diffused phase transition. High-resolution transmission electron microscopy reveals that edge dislocations play a critical role in forming locally disordered polar structures. When integrated it into a memristor device, the relaxor crystal with nanodomains demonstrates an increased number of analog resistance states and lower operating voltages compared to conventional ferroelectric crystals. This work establishes a phase engineering strategy to tailor ferroelectric behavior in low-dimensional materials, offering substantial potential for future electronic applications.

Introduction
Significant research has been focused on polarization switching in ordered ferroelectric materials, which plays a key role in advanced electronic devices, such as ferroelectric field-effect transistors, ferroelectric tunnel junctions, and ferroelectric memristors1–3. The urgent demand for high-density integration and ultralow-power consumption has prompted the exploration of ferroelectric thin film at atomic thickness limit4,5. However, sustaining ferroelectric order at reduced thickness remains a critical challenge, as it is severely perturbed by delicate energy balance that is driven by boundary conditions and depolarization field6,7. Although ferroelectrics have been achieved in thin films via epitaxial growth on lattice-matched substrates, the dependence on such substrates constrains integration into electronic devices, and thus limits potential practical applications8,9. In contrast, two-dimensional (2D) van der Waals (vdW) ferroelectric materials, featured by the absence of dangling bonds and remain the ferroelectricity down to monolayer thickness, offer a transformative platform for next-generation miniaturized electronic devices10–12.

Controlling the phase in 2D vdW materials allows probing exotic physical phenomena in condensed matter physics, particularly enabling the full exploitation of 2D ferroelectrics at the atomic thickness13–15. Considerableresearch has been focused on phase engineering for creating 2D heterostructures, superlattices and moiré systems through straightforward mechanical exfoliation and restacking, that reveals distinct physical properties and exceptional performance16,17. While this approach has led to ground-breaking discoveries in sliding ferroelectrics, topological polar textures and superconductivity, the scalability and practicality of the assembled heterostructures remain significant limitations for technological implementation18,19. Capturing the full potential of 2D materials requires robust and scalable synthesis of these atomical scale thickness materials and their heterostructures with designable spatial modulation of chemical compositions and electronic structures. The extreme aspect ratio, lack of intrinsic substrate and the atomically thin nature render them highly susceptible to defects and external perturbations (electronic, optical and mechanical)20,21. Thus, precise and efficient phase engineering is a critical challenge for unlocking the transformative potential of low dimensional materials.

Analogous to non-stoichiometric phase engineering in conventional oxides, the multiple phases in 2D materials can be finely synthesized by tailoring the composition ratio22,23. Here, we demonstrate a mixed-phase 2D system exhibiting a ferroelectric phase transition by controlling the ratio between monoclinic and trigonal phases that coexist in CuInP2(S1-xSex)6 crystals. The phase control could induce superlattice modulations and local heterogeneous polar structures, facilitating the transition from ferroelectric (FE) to relaxor ferroelectric (RFE), superparaelectric (SPE), or paraelectric (PE), which is manifested by flatten free energy profile. The synthesis and composition of CuInP2(S1-xSex)6 was characterized using energy dispersive X-ray spectroscopy (EDS), electron probe microanalyzer (EPMA), Raman spectroscopy. The local heterogeneous polar structure was investigated by X-ray diffraction (XRD), high-angle annular dark-field scanning transmission electron microscopy (HAADF-STEM), selected area electron diffraction (SAED), and second harmonic generation (SHG). Temperature-dependent ferroelectric and dielectric measurements confirmed the relaxor behavior through a hysteresis-free polarization response and frequency dispersion of dielectric permittivity. Additionally, the evolution of domain size and coercive voltage was investigated using piezoresponse force microscopy (PFM). Finally, the integration of 2D relaxor crystals into memristor devices exhibit enhanced analog resistance states and low operation voltages, demonstrating the potential for advanced electronic applications.

Results
Strategy for phase-engineering
The introduction of external elements into crystals has been an effective strategy for tuning structural phase transitions and material properties. In 2D layered transition-metal thio- and seleno-phosphates, CuInP2S6 crystallizes in ferroelectric monoclinic phase (space group Cc), whereas the CuInP2Se6 adopts a distinct ferroelectric trigonal phase (space group P31c)24,25. Inspired by this structural diversity, we hypothesized that partial substitution of S with Se in CuInP2(S1-xSex)6 crystal framework could lead to a controlled coexistence of these distinct phases and phase boundaries. Such phase coexistence results in local heterogeneous structure and polar nanoregions, which are expected to give rise to novel ferroelectric behaviors26,27 (Fig. 1a).

To systematically explore the evolution of structure and ferroelectric properties, we synthesized a series of CuInP2(S1-xSex)6 crystals using the chemical vapor transport (CVT) method (details are presented in the Experimental Section). The Se substitution is confirmed by EDS and Raman analyses, and the accurately Se content of 0%, 0.8%, 2.6%, 4.2% and 6.9% are obtained using EPMA (Supplementary Note 1, Supplementary Fig. 1 and Supplementary Table 1). Given that Se has a larger ionic radius and lower electronegativity compared to S, the substitution introduces both structural and electronic and structural perturbations to the lattice, weakening electrostatic constraints of Cu1+ cation and consequently altering the hopping motion of ions and structural stability of crystals25,28. As Se content increased from 0% to 6.9%, the gradually increased vibration peak at 346 cm−1 in Raman spectroscopy demonstrate the emergence of [P2(S1−xSex)6]4− anion (Fig. 1b). In contrast, the intensities of the Cu1+ vibration peak at 315 cm−1 is reduced with higher Se content, indicating that the mixed S-Se chalcogen environment disrupted the hopping motion of the Cu1+ cations. And the simultaneous broadening of the both 315 cm−1 and 346 cm−1 peaks demonstrates that Se substitution enhance lattice instability, which drives the formation of heterogeneous structures. These compositional modifications significantly influence ferroelectric behavior, which is revealed in phase diagrams derived from the temperature-dependent ferroelectric and dielectric permittivity analysis (Fig. 1c, Supplementary Note 2 and Supplementary Figs. 2 and 3). As the Se content increases, the ferroelectric transition from FE to RFE and subsequently to SPE phase was observed at room temperature, and the Tm (where the dielectric permittivity reaches its maximum) undergoes a significant shift toward low temperature. At low temperatures, the FE converts to antiferroelectric (AFE) phase only. Notably, these phase transitions occur at the Se content of 2.6%, which suggests that the phase coexistence is established at this composition.

Local heterogeneous polar structure
To elucidate the mechanism underlying these phase transition, the microstructure of the crystal was performed using XRD and SAED. The clear diffraction peaks at 13.4°, 27.1°, 41.3° and 56.3° corresponding to the (002), (004), (006) and (008) planes, respectively, which indicate the 2D crystal stacks along the c-axis (Fig. 2a). Notably, for x=2.6%, the local substitution of Se induced gradient strain leads to grain size refinement, interlayer dislocations, and amorphous phases, causing the broadening of the (002) diffraction peak29. As the Se content surpass 2.6%, the more uniform distribution of Se effectively releasing the accumulated strain and improving the local stacking order, resulting in the recovery of sharper (002) diffraction peak. And the diffraction peaks shift to larger angles, reflecting a reduction in interlayer distance along the c-axis. Additionally, the gradual decrease in intensity of the (004), (006) and (008) peaks with increasing Se content is related to local heterogeneous structure. These structures modifications might result in the emergence of trigonal phase with small interlayer distance, where the S atoms shift towards to the layer gap, and Cu atom move towards to the inner layer induced by steric hindrance effects.

Furthermore, based on SAED patterns along the (001) axis, the typical diffraction spots of monoclinic structure were identified at x=0%, consistent with prior reports30,31 (Fig. 2b). For x=2.6%, the superlattice diffraction spots corresponding to the monoclinic phase diffraction, along with weak trigonal phase diffraction spots demonstrate the phase coexistence structure. The corresponding atomic-scale structural evolution can be directly observed in the HAADF-STEM images along the (001) direction (Supplementary Fig. 5). The superlattice structure was confirmed through the satellite peaks in the reciprocal lattice at bc plane for x=2.6%, which was investigated by single-crystal X-ray diffraction32 (Supplementary Fig. 6). These results indicate that the integration of trigonal and monoclinic phases forms a superlattice structure within the 2D crystal, leading to a locally heterogeneous polar structure.

Furthermore, the HAADF-STEM was conducted to uncover the atomic-scale microstructure of the coexisting multiphase. During the CVT growth, the Se atoms incorporate into the lattice under near-equilibrium conditions, leading to macroscopically homogeneously distribution at S sites. The monoclinic and trigonal phases can be determined by the angles between a- and c-axis and the interlayer distance along the (010) plane (Supplementary Note 4 and Supplementary Fig. 7). For Se content below 2.6%, the substituting Se atoms are randomly and sparsely distributed within the 2D layers, and the induced point defects are insufficient to distort the lattice, maintaining a homogeneous single-phase structure. And the independent monoclinic phase was observed, characterized by the angle of 98.9° between the a and c axis and clearly vdW gap (Supplementary Figs. 5 and 8). With the Se content increase to 2.6%, local Se enrichment near defect sites induces sufficient strain to distort the lattice, leading to the formation of dislocation and trigonal phase. These localized strain fields further induce interlayer sliding within the layered crystal structure, leading to the coexistence of monoclinic and trigonal phases, as clearly observed in Fig. 2c. Edge dislocation regions, marked by orange dashed lines, indicated that substitution of S by Se, with its larger atomic radius, produced lattice stress by steric hindrance effects as illustrated in the lower right section of Fig. 2c. This can be further corroborated by the presence of dislocation cores identified through the geometric phase analysis (GPA) featured by pronounced shear stress (εxy) concentration. The strain gradient facilitated the formation of trigonal phase under compressive strain between the dislocation and monoclinic phase, resulting in nearly 90° angle between the a and c axis in unit cell. Multiphase structure that comprises dislocation regions (~5 nm), trigonal phase regions (~20 nm), and residual monoclinic phases, formed disordered polar regions (Supplementary Fig. 8), which provides detailed insights into the local polar nanoregions in the multiphase coexistence structure.

The SHG signal is a sensitive probe for detecting the local symmetry breaking in materials, with the intensity (I) is proportional to the square of the polarization (P2). For x=0%, the uniform monoclinic phase exhibit strong SHG signals due to pronounced asymmetry and large polarization ordering (Supplementary Fig. 9). As the Se content increases, the emergence of the trigonal phase disrupts polarization ordering and reduces asymmetry, which leads to the gradual decline of SHG intensity. The angle-dependent SHG measurements reveal similar sixfold rotation symmetry patterns across samples, even for the weakest intensity crystal at x=6.9%, indicating the persistence of polarization in the heterogeneous polar structure (Fig. 2d). Besides, the SHG mapping exhibits uniform polarization regions for x=0%, while the small and weak polarization region was observed at the x=2.6%, corresponding to the coexisted phases with local heterogeneous structure (Fig. 2e). The detailed analysis of polarization region areas shows that most polar regions for x=2.6% are concentrated around 0.2 μm2. These results confirm that the coexisting multiphases induce local heterogeneous structure, which decreases polar anisotropy and forms polar nanoregions.

Relaxor ferroelectric property
The coexisted multiphase with local heterogeneous polar structure is closely associated with the ferroelectric behavior. The polarization-electric field (P–E) loops and corresponding current-electric field (I–E) curves of the mixed-phase crystals CuInP2(S1−xSex)6 were measured at room temperature (Supplementary Note 5 and Supplementary Fig. 10). To avoid leakage current influence, the P–E loops was carried out at frequency 800 Hz (Fig. 3a). At x=0%, the sample exhibits typical ferroelectric feature with large maximum polarization (Pm) of 4.4 μC cm−2 and visible remanent polarization (Pr) of 3.0 μC cm−2, which leading to clearly hysteresis loops. As the Se content increases to x=2.6%, the Pr dropped sharply to 1.2 μC cm−2, while the Pm remains relatively stable at 3.7 μC cm−2. The reduced hysteresis loops indicate transition from FE to RFE state. With further increases in Se to 6.9%, both Pr and Pm decreased to 0.2 and 1.0 μC cm−2, respectively, which displays slim and nearly linear hysteresis loops characteristic of the SPE-like state. The similar polarization evolution of Pr and Pm is also observed in the exfoliated few-layer crystals, confirming that the relaxor ferroelectric behavior can be preserved in this two-dimensional crystal (Supplementary Fig. 11). Notably, the crystal of CuInP2(S1−xSex)6 possesses a ferrielectric configuration, which the Cu and In ions exhibit antiparallel displacements with unequal magnitudes, resulting in a non-zero net polarization33. Since the contribution of In ion to the polarization is negligible, the ferroelectric behavior is fundamentally governed by the displacement and positioning of Cu ions, which can occupy three types of position31,34: off-center position (Cu1), center position (Cu2); and the interlayer space (Cu3) (Fig. 3b). In the monoclinic phase, Cu ions at S plane predominantly occupy the off-center Cu1 site, which leads to a significant polarization due to large ionic displacement. Conversely, in the trigonal phase, Cu ions preferentially migrate towards the more centered Cu2 intralayer position, resulting in weaker polarization due to a much small displacement. Density functional theory (DFT) calculations also confirmed that the trigonal phase (CuInP2Se6) has weaker polarization compared to the monoclinic phase (CuInP2S6) (Fig. 3c). Moreover, the trigonal phase was found to have lower energy barriers, which facilitates polarization reversal (Fig. 3d). These results collectively confirm the FE-RFE-SPE transitions in the mixed-phase crystal with the increase of Se content.

Another feature of relaxor ferroelectric is the presence of a diffused phase transition and frequency dispersion around the temperature-dependent dielectric peaks (Tm). To characterize these properties, the dielectric permittivity of mixed-phase CuInP2(S1−xSex)6 crystal was measured across a temperature range of −150 °C to 100 °C and frequency range of 300 Hz to 3 MHz (Supplementary Note 3 and Supplementary Fig. 3). In contrast to the sharp dielectric peaks for x=0% with typical ferroelectric phase, the gradual shift of phase transition temperature at relatively small range with frequency confirms the relaxor-like ferroelectric at composition x=2.6% and 6.9% (Fig. 3e). The degree of dispersion was quantified using the relaxor diffuseness factor (γ) derived from modified Curie-Weiss law35,36: 1/ε − 1/εm = (T − Tm)γ/Cm, where the εm is the maximum dielectric constant at Tm, and Cm is the Curie-like constant, where the γ=1 corresponds to a typical ferroelectric behavior and γ=2 represents an ideal relaxor ferroelectric. The γ values increase progressively from 1.11 for x=0% to 1.45 for x=2.6%, and 1.9 for x=6.9%, demonstrating the enhanced relaxor property (Fig. 3f). And the γ of 1.45 at x=2.6% indicates the initial transition from FE toward to RFE, which corresponding to the emergence of phase coexistence structure. Notably, for x=2.6%, an abnormal behavior of the dielectric peak Tm shifts towards lower temperatures as frequency increases was observed, which attributed to the ionic conduction effects previously reported in relaxor ferroelectric materials29,37 (Supplementary Note 3 and Supplementary Fig. 4a). Besides, the gradually decrease of phase transition temperature Tm was observed as the Se content increases, which could be attributed to the local heterogeneous structure and lattice instability. At x=2.6%, the Tm is observed at 25 °C, close to room temperature (Supplementary Fig. 9b). Consequently, the maximum dielectric permittivity was observed at x=2.6%, with the subsequent decrease of dielectric permittivity for the x exceeding 2.6% at room temperature (Supplementary Fig. 9c). Moreover, the thermal stability of dielectric permittivity was substantially improved with the greater Se content. For the x=6.9%, the highly stabilized permittivity of approximately 60 was maintained across a wide temperature range from −100 °C to 50 °C, in contrast to x=0%, where stability was limited to a narrower range (−100 °C to 20 °C) with a permittivity of 14. This enhancement is attributed to the pronounced relaxor behavior, which enables the polarization to gradually vanish even above the Tm. And the increase of dielectric permittivity is associated with the local polar regions, which allows larger polarization response at lower electric field. Additionally, as the temperature approached the relaxor freezing temperature (Tf), a phase transition from AFE to RFE was observed, as depicted in Fig. 3g and Supplementary Fig. 12, further corroborating the phase evolution shown in Fig. 1c. These results demonstrate the appearance of relaxor ferroelectric at the phase boundary with monoclinic and trigonal phase.

Electrical properties
The appearance of local heterogeneous polar structure typically causes the domains size refinement within the phase coexistence. As validated by PFM, the clear reduce of the domain size is observed in few-layer crystals as the Se content increase (Fig. 4a, Supplementary Note 6 and Supplementary Figs. 13 and 14), with domains nearly disappear at the x=4.2% and 6.9%. This further demonstrate that the relaxor ferroelectric behavior could be remained in the low dimensional crystal. According to the Landau phenomenological theory, the competition between the refined domains energy and interfacial energy lowers the energy barrier for polar rotation, and makes the flatten free energy profile27. This was further corroborated by PFM phase hysteretic loop, which showed decreasing switching voltages for domain reversal with increasing Se content (Fig. 4b). Therefore, the formation of nanodomains in relaxor ferroelectric enables multiple resistance state modulation under lower operating voltages, making it well-suited for the memristor devices rely on repeated write-read operations. To further investigate the resistance behavior during the domain-switching process, the crystal was exfoliated and integrated into two-terminal memristor device with a vertical structure (Au/Ti/CuInP2(S1−xSex)6/Au), which could mimic artificial synaptic functions (Supplementary Note 7). From the sweeping voltage from −3 to 3 V at a 0.03 V step, the clear hysteretic current-voltage (I–V) curves were observed between the high resistance state (HRS) and low resistance state (LRS) for 0%, 2.6% and 6.9% (Fig. 4c and Supplementary Fig. 15). The maximized switching ratio (ILRS/IHRS) was achieved at write voltage of 1.62 V for x=2.6%, due to the heterogenous polar structure enhanced polarization rotation in the relaxor ferroelectric structure (Fig. 4d). Furthermore, the resistance switching dynamics, linked to polarization rotation, was analyzed using pulse sequences voltage (Supplementary Fig. 16). The upward and downward polarization were defined as ON and OFF resistance state, respectively. The device was first set to the OFF state by pulse voltage Vset = –3 V (duration t = 50 ms), followed by the application of pulse voltages Vwrite and Vread. Over the 40 consecutive pulses, the resistance exhibited step-wise increases, mimicking the long-term potentiation (LTP) behavior of artificial synapses (Fig. 4e). Notably, the multi-state resistance of RFE at 2.6% keeps increasing depending on the pulse number, while the FE and SPE with x=0% and 6.9% reaches saturation rapidly after a few pulses. This indicated the RFE exhibits many analog resistance switching states in the memristor, which can be ascribed to multiple polarization states as nanodomains gradually switching. Additionally, these refined nanodomains effectively reduce the coercive field in RFE crystals, making low operation voltages for memristor device. At x=2.6% and x=6.9%, the required lowest Vwrite = 3 V and Vread = 1.5 V were lower than x=0% with Vwrite = 4 V and Vread = 2.5 V. These low operation voltages making relaxor ferroelectric crystal candidates for highly energy-efficient memristor applications.

Methods
Growth of 2D crystal
The single crystal of CuInP2(S1-xSex)6 were growth by the chemical CVT method with iodine as the transport agent38. High pure precursor materials, including Cu (99.99%, Sigma-Aldrich), In (99.99%, Sigma-Aldrich), P (99.999%, Sigma-Aldrich), S (99.98%, Alfa Aesar) and Se (99.99%, Sigma-Aldrich) were mixed in quartz tube with a nominal stoichiometric proportion of 1:1:2:6:0, 1:1:2:5.982:0.018, 1:1:2:5.916:0.084, 1:1:2:5.85:0.15, 1:1:2:5.748:0.252, corresponding to the Se content 0%, 1%, 3%, 5% and 8%, respectively. This slightly deviation is ascribed to the variations in vapor transport efficiency and minor weighing inaccuracies of the elemental precursors during the CVT growth process. Then they were sealed under a vacuum of 10−4 mbar through the high-temperature torch. The crystals were growth in a dual-temperature zone tubular furnace, with source zone maintained at 750 °C and the growth zone at 700 °C, over a period of 7 days. The gradually heating and cooling temperature rate of 2 °C min−1 was maintained throughout the process. Upon completion of the growth, plate-like crystals with lateral dimensions of several millimeters and thicknesses of several micrometers were obtained in the crystallization zone.

Structure characterization
Detailed crystal structure analysis was conducted using Single-crystal X-ray diffraction with Bruker D8 Venture diffractometer to determine the crystal structure, with data processing, cell refinement, and calculations conducted via APEX4, SAINT, and SADABS software. And the crystals were mounted on silicon during XRD measurements. Elements concentration was characterized using a scanning electron microscope (Nova SEM) equipped with an EDS and EPMA (JEOL JXA-iHP200F). Raman spectroscopy was acquired using a WITec alpha300 apyron confocal system with a 532 nm excitation laser. Cross-section TEM lamellae were prepared using a Thermo Fisher Helios G4 system. Atomic-scale conformations were captured with HAADF-STEM using an FEI Titan ST transmission electron microscope.

SHG measurement
Polarization-dependent SHG measurements and SHG mapping were conducted using a WITec alpha300 Raman system. A 1064 nm excitation laser beam was focused on the objective focal plane, with its polarization controlled by a half-wave plate. Angle-dependent SHG measurements were performed in 2° increments using a motorized rotational stage, and the SHG intensity at each angle was averaged over three counts with an integration time of 0.5 seconds per step. For SHG mapping, a selected region of 10 μm × 10 μm was scanned under consistent incident light intensity, with a scan rate of 50 ms per pixel.

PFM measurement
Ferroelectric domain characterization was conducted using piezoresponse force microscopy (PFM) on an Asylum Research MFP-3D system, operating in dual AC resonance tracking mode (DART-PFM) at room temperature. Conductive Pt-coated Si tips with a spring constant of 3 Nm−1 and a resonance frequency of 75 kHz were used for the measurements. Domain imaging was performed with an AC voltage of 1 V applied at the tip, operating at a tip-sample contact resonance frequency of 315 kHz. The phase switching voltage was measured by applying the sweep voltage from −9 to 9 V with 1 Hz and drive voltage of 1 V at the tips, the phase curves can be obtained at the off state.

Device fabrication and measurement
The 2D crystal flakes were obtained by mechanically exfoliating them from the bulk crystal and were identified as high-quality using optical microscopy. The patterned bottom electrode was fabricated by depositing 10 nm Ti and 30 nm Au onto a SiO2/Si substrate using a laser writer (DWL66+) and e-beam evaporation. The exfoliated 2D crystal flakes were then transferred onto the bottom electrode using the dry transfer method with a polymer stamp. The top electrode was defined using a laser writer, followed by deposition of 10 nm Ti and 80 nm Au through e-beam evaporation. Electrical measurements of the device were performed using a Keithley 4200 semiconductor characterization system at room temperature in an ambient environment.

Ferroelectric and dielectric measurement
Before the measurement carried out, the circle Au electrode (diameter 1 mm) was deposited on both sides of the crystal at thickness of 100 nm. The ferroelectric property was investigated by the ferroelectric analyzer with a Sawyer-Tower circuit (axiACCT, TF Analyzer 2000E). The P–E loops were performed by the bipolar electric fields of the triangular wave in DHM mode, and the frequency from 100 Hz to 1 kHz. The temperature dependent P–E loops were obtained from –120 to 10 °C at increase step of 10 °C. The dielectric permittivity was acquired over a broad frequency using the Alpha-A High Performance Frequency Analyzer by applying an AC voltage of 0.1 V. With the temperature increase from –150 to 100 °C, the dielectric spectroscope from 300 Hz to 3 MHz was measured at the temperature increase step of 2 °C.

DFT calculations
First-principles calculations were performed using the Vienna Abinitio Simulation Package (VASP) with the projector-augmented wave (PAW) method and the Perdew-Burke-Ernzerhof (PBE) exchange-correlation functional39,40. A plane-wave kinetic energy cut-off of 450 eV was chosen. For structural relaxations, a Gamma-centered k-points grid with ka > 40 Å was used to sample the Brillouin zone, where k is the number of mesh points and a represents the lattice constants. A similar k-point mesh was used for all static runs. To prevent any spurious interaction between periodically repeated layers, a large vacuum space of at least 20 Å in the direction of c was applied. Structures were relaxed until the energy differences converged below 10−6 eV and forces within 1 × 10−2 eV/Å. van der Waals corrections using the optB86b exchange-correlation functional were applied to precisely account for interlayer interactions, particularly for materials like CuInP2S641,42. Solid-state nudged elastic band method (SSNEB) was employed to investigate the polarization reverse barrier43. The spontaneous polarizations of the ferroelectric materials were computed using the Berry phase method, which incorporates contributions from both electronic and ionic parts44.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "The authors present their work on a phase engineering approach of van der waals crystals which allowed them to achieve and demonstrate the tailoring of the ferroic behavior of CuInP2S(1-x)Se(x)6 from ferroelectricity to relaxor ferroelectricity to superparaelectricity. The authors postulate and demonstrate that substitution of S atoms with Se will lead to a mixed phase system of monoclinic and trigonal phases which can reveal new and interesting physical properties. They were able to show a mixed phase diagram with support from TEM imaging to show the presence of possible polar nanoregions due to the phase mixture, which they attribute as the main contributing factor to the observation of the relaxor phase. Although the results seem to have some novelty, a careful evaluation of this manuscript raises several questions that must be first answered before it can be considered for publication in a high impact factor journal such as Nature Communications. The following lists my main concerns:\n1. Novelty of the results: The authors claim to have demonstrated relaxor ferroelectrics in low-dimensional materials for the first time (line 20). There seems to be a convolution between the exact meaning of the terms van der waals crystals (which are typically bulk crystal) and low-dimensional systems such as monolayers or few layer thick materials. This is important because as Song et al (Phys. Rev. B 96, 235420, 2017) have shown, the ground state of monolayer CIPSe can be very different from the bulk CIPSe crystals. In addition, the phase engineering of the multiple phases has been demonstrated earlier in bulk crystals of CuInP2SxSe(1-x), and phase diagrams of the coexistence of mixed phases and mixed ferroic behavior - relaxor as well as dipole glassy behavior has been demonstrated in the work by Macutkevic et al. (Phys. Rev. B 78, 064101, 2008). If the current results were demonstrated in a monolayer or few layer thick geometry, the results could have been considered as novel. Can the authors demonstrate such relaxor behavior in a monolayer or few layers thick systems and will the phase diagram remain the same in such systems?\n2. Technical questions:\ni. The phase diagram in Fig. 1c was said in the manuscript to be constructed from a combination of the temperature dependent P-E loops and dielectric permittivity analysis in supplementary Fig. S2 and note S2. However, figure S2 contains only the P-E loop information wrt temperature. This figure raises several questions. First of all, the x-axis of Fig. S2 is highly confusing. Second, these P-E loops were taken up to 10C but in the phase-diagram of Fig. 1c, there are data points up to Room temperature (25C). Can the authors show the RT data as well? if the loops are leaky, maybe a higher frequency might be considered. Third, it is not immediately clear from the broadening of the loops with temp and corresponding increase in remanent polarization for x=1.4% how they arrived at the conclusion that the loops resemble relaxor-like? Fourth, what is the origin of the difference in the shape of the low temperature AFE loops for x=1.4% and x=2.5% with the x=2.5% loops resembling more of the typical AFE-like loops?\nDid the authors mean to also include the information from supplementary note 5 and corresponding Fig. S8 to help explain Fig. 1c? The phase diagram should be explained better - it is confusing currently as to how it was arrived at.\nii. In the supplementary figure 8, the peaks of the real dielectric constant seem to shift towards a lower temperature with an increase in the frequency, but the typical dielectric behavior is opposite - the peak shifts towards higher temperatures with increasing frequency. The authors refer to reference 8 in the supplementary to explain this discrepancy. But, from as far as I can tell Fig. 3 in Ref 8 has the opposite trend to what the authors show here. Can the authors explain this discrepancy?\niii. Is it possible to estimate the phase fraction of the trigonal phase with increase in the Se doping concentration?\niv. The authors highlight the importance of the strain gradients under Se substitution towards the formation of the trigonal phase. My question in this regard is related to point 1 above. How feasible will it be to generate such strains during the growth of a just a few monolayer thick films due to the differences in the mechanical boundary conditions between crystal growth vs film growth?\nv. What is the physical and technological significance of the relaxor diffuseness factor with Se doping of x=1.4% being 1.45 mean? At x=4.2%, and the diffuseness factor of 1.9 seems to suggest a much stronger relaxor-like behavior, but this phase is stabilized only at low temperatures but for most practical technological applications, the temperature is usually RT or above.\nvi. Can the authors elaborate more on why they considered relaxors for potential memristor applications? Relaxors are known to have lower remanent polarization which would mean that they are more susceptible to retention loss which would seem detrimental for any type of memory applications. It might have been more useful to quantify the electrostrictive coefficients or measure the strain in such relaxor systems.\n3. Other comments:\ni. Please use a consistent system of labeling the supplementary figures either calling them Supplementary Figure x or Figure Sx.\nii. Fig. 3g: the figure says the transistion is from AFE to FE but the label says it shows transition from AFE to RFE?\niii. Fig.2c: labels of atoms are missing", "reviewer_2": "The manuscript by Yang et al. is entitled “Phase engineering of relaxor ferroelectricity in van der Waals crystal” (id: NCOMMS-25-70669-T), which realizes the relaxor ferroelectricity in low-dimensional materials by controlling phases. The phase engineering is realized by partial substitution of S with Se in CuInP2(S1-xSex)6 and establishes competing coexistence of monoclinic and trigonal phases. This approach reveals a sequential transformation of FE to RFE and SPE phases at room temperature in the 2D materials. Through experimental techniques like XRD, SAED, TEM images, SHG signal, and PFM, they tried to establish their goal. The data quality seems sufficient to support the main conclusions. They use DFT to compare the energy barrier of phase transition and polarization of the two pristine monoclinic and trigonal phases.\nMost of the manuscript is well written; however, it would benefit from a careful revision to eliminate repetitions, correct typos, and fix mislabeled references. Before recommending publication in Nature Communications, the authors should consider the following points.\ni) The authors claim that they have realized the relaxor ferroelectricity in 2D materials for the first time. However, there is a report that previously indicated the realization of relaxor ferroelectricity in Se-doped 2D PbTe by Gong et al. (Strain-induced two-dimensional relaxor ferroelectrics in Se-doped PbTe; DOI: https://doi.org/10.1103/PhysRevB.109.054117). It is not clear why the authors insist on this statement. Either the authors should be clearer with their statement or remove such a statement.\nii) As mentioned, the XRD diffraction peaks at 13.4° correspond to the (002) indicates the 2D crystal stacks along the c-axis. According to the Authors, broadening the (002) diffraction\npeak suggests that atomic substitution induces grain size refinement, interlayer dislocations, and\namorphous phases for Se content of x = 1.4%. However, when the concentration of Se further increased to x = 2.5% and x = 4.2%, the (002) peak again became sharp in nature. The reason for this anomaly is not discussed in the manuscript. Authors are advised to discuss this to clear the anomaly.\niii) Since these are the multilayers of CuInP2(S1-xSex)6. There is a high possibility of having different concentrations of Se substitution in the different layers, leading to differences in the lattice parameter. Also, as the coupling between the layers is expected to be weak, a moiré pattern should emerge; did the authors check for it with STM? If they have observed that, then what is the role of that varying interlayer coupling for this phase transformation?\niv) Using DFT, the Authors tried to support the easy transition along FE-RFE-SPE. However, the results of Fig. 4(c) and 4(d) are for pristine CuInP2Se6 and CuInP2S6. Furthermore, this result is derived by connecting the P_down to the P_up structure. The easy transition along FE-RFE-SPE, using these pristine calculations, does not make strong evidence for the experimental support. The calculations using CuInP2(S1-xSex)6 structure will give more realistic evidence.\nv) The Authors mentioned that there is a decrease in switching voltages for domain reversal with increasing Se content (Fig. 4b). This claim was valid and consistent from x = 0% and x = 2.5%. However, for x = 4.2%, the switching voltage is higher than x = 1.4% and 2.5% on the positive side. Please clear this anomalous observation.\nMinor comment:\ni) The texts in SI from line numbers 66 to 106 have an issue. Also, after line number 106, line number 265 has started.", "reviewer_3": "In this manuscript, the authors addressed an interesting work on phase engineering of relaxor ferroelectricity in van der Waals crystals. Authors demonstrated a mixed-phase 2D CuInP2(S1-xSex)6 single crystals where tuning the ratio of monoclinic and trigonal phases induced superlattice structures and heterogeneous local polar states, enabling transitions from ferroelectric to relaxor, superparaelectric, or paraelectric behavior through a flattened free-energy profile. Structural, compositional, electrical, and atomic scale measurements were performed using various characterization techniques. Temperature-dependent dielectric and ferroelectric measurements confirmed relaxor behavior via hysteresis-free polarization and frequency dispersion of permittivity, while PFM revealed domain size and coercive voltage evolution. Integration into memristor devices further showed enhanced analog resistance states and low operating voltages, highlighting potential for electronic applications. Although the manuscript was interesting. However, I regret that the explanation and mechanism of the manuscript were not convincing to support many of their major claims. The mechanism of the study remained confusing due to incomplete quantitative evidence and limited analysis. In its current form, the manuscript cannot be considered for publication until all the raised concerns were properly addressed.\nMajor comments\n1) A central concern was the issue of composition control and quantification. The reason for selecting a non-uniform substitution trend of the Se ratio into the host lattice was unclear. Typically, to establish compositional trends and phase boundaries, systematic increment doping strategies should be employed. The current approach raised questions about how the authors controlled and verified the exact Se content, particularly at the extremely low substitution (around 1.4%) levels reported. Without clear justification and compositional control, the scientific basis of the correlation between composition and physical behavior was unclear.\n2) The authors synthesized a series of CuInP2(S1-xSex)6 crystals with the increasing Se content of 0%, 0.3%, 1.4%, 2.5%, and 4.2%. Based on the formula, the Se should be homogenously substituted rather than the formation of a composite single crystal, where a mixed region of monoclinic and trigonal would exist with shared boundaries, where two regions apply interfacial strain. In composites, lattice mismatch at interfaces could generate strain. The authors should discuss why uniform substitution was pursued over this alternative and whether the substitution was homogeneous or whether the authors synthesized the composite single crystals.\n3) In Supplementary Table 1, the atomic percentage values of Cu, S, and Se were inconsistent, which challenges the compositional analysis. EDS alone was insufficient for quantitative determination, especially at such low substitution levels. The authors should validate their results with more precise techniques, such as the electron probe microanalyzer (EPMA) technique [Phys. Rev. Materials 3, 084405 (2019)].\n4) Phase identification was another weak point of the manuscript. The authors showed the coexistence of monoclinic and trigonal phases, and the proposed superlattice formation, qualitatively. But there was no quantitative presentation, like reporting of lattice parameters or phase fractions.\n5) Based on the STEM results, the authors attributed the formation of edge dislocations to the Se substitution, but the doping levels were too low to justify such noticeable defect generation. Alternative explanations, such as strain gradients at grain boundaries, temperatures, or intrinsic growth-related stresses, must be considered.\n6) Moreover, the authors correlated the relaxor behavior with the formation of the edge dislocations, which again raised concern about the clear mechanism of the current work. A systematic investigation was required to establish whether these dislocations were truly induced by substitution or were inherent to the growth process using quantitative defect density analysis.\n7) The SHG results presented for the pure composition show dark regions that are not explained. Since SHG signals typically appeared bright in non-centrosymmetric regions where inversion symmetry was broken [Adv. Mater. 35, 2210562 (2023)], the presence of dark regions in nominally ferroelectric domains was relatively confusing. The authors need to discuss whether these dark regions arose from structural inhomogeneity, phase coexistence, or during the experimental process.\n8) Peak splitting in current-electric field loops provides direct evidence of antiferroelectric switching. As the switching current measurements are accessible to the authors. They should include these results to evidence the presence of antiferroelectric-like behavior by showing the two switching peak characteristics.\n9) The P-E hysteresis loops exhibited noticeable imprint behavior, but the origin was not discussed. Given the presence of edge dislocations, it was plausible that internal fields from defect dipoles or strain gradients contribute to double-like hysteresis loops, where pinched hysteresis loops and linked double switching current behaviors were observed with the formation of dislocations [ACS Appl. Electron. Mater. 6, 8328-8338 (2024)]. The authors should analyze and provide further discussion on whether the imprint was correlated with dislocation/defect densities.\n10) In Raman spectroscopy (Fig. 1b), the authors claimed that the intensities of the Cu1+ vibration mode at 315 cm-1 weaken with higher Se content, indicating that the mixed S-Se chalcogen environment disrupted the hopping motion of the Cu1+ cations. However, on a close look, it seems like the S-P-S and Cu+1 peaks were more widened with the increasing Se content, with Cu+1 nearly comparable in intensity to the intensity of the peak at 346 cm-1.\n11) In Figure S2, the labeling of the horizontal axis (electric field vs. voltage) was unclear.\n12) For fair presentation of the results, the authors should show the full XRD results in wider range.", "reviewer_4": "Tao Yang et al synthesized layered van der Waals CuInP2(S1-xSex)6 crystals with controllable coexistence of the monoclinic and trigonal phases, which in turn leads to the appearance of relaxor-like ferroelectricity in the phase coexistence region. This, indeed, can transform the ferrielectric long-range order of CuInP2(S1-xSex)6 into relaxor-type ferroelectricity with superparaelectric features. As evidence of relaxor-type ferroelectric properties, the hysteresis-less polarization switching and diffused phase transition were demonstrated. Also, HR-TEM confirms the principal role of edge dislocations play in formation and stabilization of locally disordered polar structures. In overall, I support authors’ claim that a phase engineering of layered van der Waals ferrielectrics, opening their great potential for advanced nanoelectronics and related applications.\nI recommend the work for publication after a major revision. I regard that the formulation of the key findings made in this work and corresponding terminology should be reconsidered and significantly softened due to the following argumentation\n1) To the best of my understanding and knowledge CuInP2(S1-xSex)6 are ferrielectrics (FI), not ferroelectrics (FE), which is confirmed by extensive literature. See experimental and theoretical grounds in the works:\nY. M. Vysochanskii, V. A. Stephanovich, A. A. Molnar, V. B. Cajipe, and X. Bourdon, Raman spectroscopy study of the ferrielectric-paraelectric transition in layered CuInP2S6. Phys. Rev. B 58, 9119 (1998),\nhttps://doi.org/10.1103/PhysRevB.58.9119\nN. Sivadas, P. Doak, and P. Ganesh. Anharmonic stabilization of ferrielectricity in CuInP2Se6. Phys. Rev. Res. 4, 013094 (2022), https://doi.org/10.1103/PhysRevResearch.4.013094\nA. N. Morozovska, E. A. Eliseev, A. Ghosh, M. E. Yelisieiev, Y. M. Vysochanskii, and S. V. Kalinin. Anomalous Polarization Reversal in Strained Thin Films of CuInP2S6, Phys. Rev. B 108, 054107 (2023)\nhttps://link.aps.org/doi/10.1103/PhysRevB.108.054107\nThe AFE-FI-FE transitions of hysteresis loops are very clearly seen from Supplementary Fig. 2, that in my opinion, can be moved to the main text. Also, the figure shows a very significant contribution from the electric conductivity (losses) to the hysteresis loops.\n2) The ferrielectric 3-well or 4-well free energy potential can be transformed into the relaxor-like relief, but the relative dielectric permittivity (shown in Fig. 3e) is not very high (or “giant”) as it should be in true relaxors. The diffuse maximum of permittivity does not exceed the values of 200-350, being much smaller than 104-106 (as in ferroelectric relaxors). While the frequency dispersion of the permittivity is present in the frequency range 300 Hz – 3 MHz, the dispersion is much weaker than in ferroelectric relaxors. As a matter of fact, the ferrielectric with a diffuse phase transition and relaxor-like properties is obtained.\n3) I failed to understand why the author called the behavior shown in the bottom of Fig.3e as superparaelectric (SPE). Here the relative dielectric permittivity is relatively low; it increases above 40-50 oC. The increase of dielectric permittivity above 40-50 oC can be related to the thermally-activated conductivity (Maxwell-Wagner type effects). This is clearly seen from Supplementary Fig. 8 and Supplementary Fig. 9, where the imaginary part of dielectric permittivity and tangent of dielectric losses increases strongly above 40-50 oC. True relaxors and superparaelectrics (if any exist) should be dielectrics in the actual range of these properties’ observation. Also, the polarization in saturation should be very high in a relaxor that response to the external electric field as a super-dipole with a large dipole moment. However, the polarization shown in Supplementary Fig. 2 for x=4.2 % is relatively low (below 2 µC/cm2), Thus, authors rather observe some SPE-like features, but not the SPE state.\n4) The work is very well illustrated, but at the same time, the important detail shown in the figures are not commented neither in the main text nor in the Supplement."}, "author_response": "Response to Reviewers\n(Manuscript No: NCOMMS-25-70669-T)\nReviewer #1 (Remarks to the Author):\nThe authors present their work on a phase engineering approach of van der waals crystals which\nallowed them to achieve and demonstrate the tailoring of the ferroic behavior of CuInP (S\n2 1-\nSe ) from ferroelectricity to relaxor ferroelectricity to superparaelectricity. The authors\nx x 6\npostulate and demonstrate that substitution of S atoms with Se will lead to a mixed phase system\nof monoclinic and trigonal phases which can reveal new and interesting physical properties.\nThey were able to show a mixed phase diagram with support from TEM imaging to show the\npresence of possible polar nanoregions due to the phase mixture, which they attribute as the\nmain contributing factor to the observation of the relaxor phase. Although the results seem to\nhave some novelty, a careful evaluation of this manuscript raises several questions that must be\nfirst answered before it can be considered for publication in a high impact factor journal such\nas Nature Communications. The following lists my main concerns:\nResponse:\nWe are really grateful for your positive and professional evaluation of our work and\nrecognition of its novelty, particularly noting that “the results seem to have some novelty” and\n“which can reveal new and interesting physical properties.”. We have carefully addressed all\nthe perceptive comments and carried out additional experiments, and thoroughly revised the\nmanuscript accordingly. We are confident that these efforts have significantly strengthened the\npaper and that the revised version meets the high standards of this journal.\n1. Novelty of the results: The authors claim to have demonstrated relaxor ferroelectrics in low-\ndimensional materials for the first time (line 20). There seems to be a convolution between the\nexact meaning of the terms van der waals crystals (which are typically bulk crystal) and low-\ndimensional systems such as monolayers or few layer thick materials. This is important because\nas Song et al (Phys. Rev. B 96, 235420, 2017) have shown, the ground state of monolayer CIPSe\ncan be very different from the bulk CIPSe crystals. In addition, the phase engineering of the\nmultiple phases has been demonstrated earlier in bulk crystals of CuInP (S Se ) , and phase\n2 1-x x 6\ndiagrams of the coexistence of mixed phases and mixed ferroic behavior - relaxor as well as\ndipole glassy behavior has been demonstrated in the work by Macutkevic et al. (Phys. Rev. B\n78, 064101, 2008). If the current results were demonstrated in a monolayer or few layer thick\ngeometry, the results could have been considered as novel. Can the authors demonstrate such\nrelaxor behavior in a monolayer or few layers thick systems and will the phase diagram remain\nthe same in such systems?\nWe sincerely appreciate your professional and insightful comments on the relaxor behavior\nin the few layers. Following your valuable suggestion, we have confirmed the relaxor\nbehavior for few-layers crystals through the P–E measurement. The few layers crystal\ndevice was prepared by mechanical exfoliated and integrated into the two-terminal sandwich\nstructure. As illustrated in Fig. R1, the obvious hysteresis in the P–E loop at Se content of 0%\nindicates the typical ferroelectric behavior of CuInP S . When the Se content increases to 2.6%,\n2 6\nboth the coercive field and remnant polarization reduces, particularly the shape of the P–E loop\nsuggests a transition from ferroelectric to relaxor-ferroelectric behavior (Fig. R1). With a further\nincrease in Se content to 6.9, the nearly linear P–E loop reveals a transition from relaxor-\nferroelectric to superparaelectric (Fig. R1). These observations confirm the relaxor ferroelectric\nbehavior can still persist in few layers crystal. When the dimension reduces from bulk to few-\nlayers, the hysteresis loop changes from typical rectangular shape to tilted form, mainly due to\nthe enhanced depolarization. On the other hand, the PFM measurements of few layers crystal\nreveal progressively reduced ferroelectric domain sizes, indicating the formation of local polar\nnanoregions that are typical of relaxor ferroelectrics (Fig. R2). These results provide further\nevidence that the relaxor ferroelectricity could persists in few layers system.\nFig. R1 (Added as Supplementary Fig. 11 in Supporting Information) P–E loops of few\nlayers crystal with different Se contents with 0%, 2.6%, and 6.9%. Measurements were\nperformed at frequency of 300 Hz and room temperature. The insets show the device structure\nused for measurement of the few-layer crystal with thickness around 50 nm. Scale bar, 20 μm.\ncontent of 0% (ⅰ), 0.8% (ⅱ), 2.6% (ⅲ), 4.2% (ⅳ) and 6.9% (ⅴ), scale bar, 20 μm. b, PFM\ntopography, phase and amplitude images of few layers crystal with the different Se content at\n0% (ⅰ), 0.8% (ⅱ), 2.6% (ⅲ), 4.2% (ⅳ) and 6.9% (ⅴ), scale bar, 1 μm.\nTo address this comment, we have updated Supplementary Fig. 11 to include the P–E\nloops of the few-layer crystals, demonstrating that the relaxor ferroelectric behavior could be\nmaintained in the low dimensional crystal. In addition, we add “The similar polarization\nevolution of P and P is also observed in the exfoliated few-layer crystals, confirming that the\nr m\nrelaxor ferroelectric behavior can be preserved in this very thin crystal plate (Supplementary\nFig. 11).” on page 5 of the revised manuscript. And we add “As validated by PFM, the clear\nreduce of the domain size is observed in few-layer crystals as the Se content increase (Fig. 4a\nand Supplementary Note 6, and Supplementary Fig. 13 and Fig. 14), with domains nearly\ndisappear at the x = 4.2% and 6.9%. This further demonstrate that the relaxor ferroelectric\nbehavior could be remained in the two-dimensional crystal.” on page 7 of the revised\nmanuscript.\n2. Technical questions:\ni. The phase diagram in Fig. 1c was said in the manuscript to be constructed from a combination\nof the temperature dependent P-E loops and dielectric permittivity analysis in supplementary\nFig. S2 and note S2. However, figure S2 contains only the P-E loop information wrt temperature.\nThis figure raises several questions. First of all, the x-axis of Fig. S2 is highly confusing. Second,\n3\nthese P-E loops were taken up to 10 °C but in the phase-diagram of Fig. 1c, there are data points\nup to Room temperature (25 °C). Can the authors show the RT data as well? if the loops are\nleaky, maybe a higher frequency might be considered. Third, it is not immediately clear from\nthe broadening of the loops with temp and corresponding increase in remanent polarization for\nx=1.4% how they arrived at the conclusion that the loops resemble relaxor-like? Fourth, what\nis the origin of the difference in the shape of the low temperature AFE loops for x=1.4% and\nx=2.5% with the x=2.5% loops resembling more of the typical AFE-like loops?\nDid the authors mean to also include the information from supplementary note 5 and\ncorresponding Fig. S8 to help explain Fig. 1c? The phase diagram should be explained better -\nit is confusing currently as to how it was arrived at.\nResponse:\nWe truly appreciate your careful reading and fully understand your concerns in Fig. 1c.\nThe composition-dependent phase diagram was primarily constructed from the temperature-\ndependent P–E hysteresis loops and dielectric permittivity for different content of Se. Therefore,\nthe phase diagram in Fig. 1c is mainly derived from Supplementary Fig. 2 and Fig. 3,\nSupplementary Note 2 and Note 3 in the revised manuscript.\nRegarding your first question, the Supplementary Fig. 2 illustrates the evolution of\nferroelectric behavior in crystals with different Se contents as the temperature increases from\n−120 °C to 10 °C. All the corresponding P–E hysteresis loops were measured under the same\nelectric field of 20 kV mm⁻1 and frequency of 100 Hz. To better clarify this point, we have\nreplotted the Figure with temperature as the x-axis, as shown in the update Fig. R3, where\nthe P–E loops were recorded every 10 °C.\nRegarding your second question, the room-temperature ferroelectric hysteresis loop\nhas been updated in Fig. R3. With increasing temperature, the broadening of the hysteresis\nloops was caused by the enhanced leakage current, therefore, we increased the measurement\nfrequency to 800 Hz to minimize this effect. In addition, the room-temperature P–E loops have\nbeen presented in the updated Fig. S4, which was measured under an electric field of 20 kV\nmm⁻1 at a frequency of 800 Hz.\nRegarding your third question, we demonstrate the relaxor ferroelectric behavior\nthrough reduced hysteresis P–E loops and diffuse phase transition. To minimize the\ninfluence of leakage current, a high-frequency P–E hysteresis loop was measured in Fig. R3.\nAt room temperature, as the Se content increases from 0% to 2.6%, the markedly reduced\nhysteresis and remnant polarization clearly indicate the transition from a ferroelectric to a\nrelaxor ferroelectric. Furthermore, in contrast to the sharp dielectric peak observed in\nferroelectrics with 0% Se content, the broadened dielectric peak at Se content of 2.6%\n(corresponding to 1.4% Se content in previous manuscript) indicates a relaxor ferroelectric\n4\nbehavior (Fig. R4). This is further evidenced by the diffuse phase transition, that the phase\ntransition temperature the gradually shifts with the frequency increase from 300 Hz to 3 MHz.\nThese findings provide compelling evidence of relaxor ferroelectric behavior in the crystal with\n2.6% Se content.\nRegarding your fourth question, the crystal of 2.6% Se exhibits AFE-like loop with\nnoticeable imprint behavior, which can be ascribed to the edge dislocation, as reported\npreviously1,2. With the Se substitution introduces abundant dislocations, the observed pinched\nhysteresis loops primarily arise from the antiferroelectric phase and dislocation-related pinning\neffect. With the further increase of Se content, the more uniformly distributed Se effectively\nreleases the accumulated strain and reduces dislocation density, while simultaneously inducing\nnonpolar regions and refining local polar nanoregions. As a result, the crystal with 4.2% Se\n(corresponding to 2.5% Se content in previous manuscript) presents a more typical\nantiferroelectric-like hysteresis loop. This domain refinement could be verified from the PFM\nimages (Fig. R5). When the Se content exceeds 2.6%, the disappearance of observed domains\nindicates the dominance of polar nanoregions, further confirming the composition-driven\nferroelectric transition.\nFig. R3 (Updated in Supplementary Fig. 2 in Supporting Information) Evolution of\nferroelectric behavior with temperature in the range of –120 to 25 °C. P–E loops were measured\nfrom –120 to 10 °C at 10 °C intervals under the same electric field of 20 kV mm–1 and\nfrequency of 100 Hz. P–E loops were measured at 25 °C under the electric field of 20 kV mm–1 and\nfrequency of 800 Hz.\n−80 °C to 75 °C for the Se content of 0%, 2.6% and 6.9%.\nFig. R5 (In Fig. 4a in Manuscript) PFM amplitude image of domain structure evolution with\nSe content increasing from 0% to 6.9%, scale bar, 1 μm.\nReferences:\n1. Sheeraz, M. et al. Defect Engineering of Ferroelectric Hysteresis in Lead-Free\nBi (Na,K) TiO Thin Films. ACS Appl. Electron. Mater. 6, 8328-8338 (2024).\n1/2 1/2 3\n2. Lee, D. et al. Active control of ferroelectric switching using defect-dipole engineering. Adv.\nMater. 24, 6490-6495 (2012).\nto temperature and adding the room-temperature hysteresis loops with higher frequency. In\naddition, we added Supplementary Fig. 3 to better illustrate the variation of the phase transition\ntemperature, and providing a clearer understanding of the origin of the phase transitions in the\nphase diagram of Fig. 1c.\nii. In the supplementary figure 8, the peaks of the real dielectric constant seem to shift towards\na lower temperature with an increase in the frequency, but the typical dielectric behavior is\nopposite - the peak shifts towards higher temperatures with increasing frequency. The authors\nrefer to reference 8 in the supplementary to explain this discrepancy. But, from as far as I can\ntell Fig. 3 in Ref 8 has the opposite trend to what the authors show here. Can the authors explain\nthis discrepancy?\nResponse:\nThank you for pointing out this important issue. Indeed, in typical relaxor ferroelectrics,\nthe dielectric peak shifts toward high temperature with the frequency increase. This due to the\ndielectric response follows the Debye-type relaxation model, in which the dielectric constant\nreaches its maximum at ωτ = 1 (where ω is the angular frequency and τ is the relaxation time).\nAs the frequency increases, the dielectric peak shifts towards smaller τ values. According to the\nArrhenius relationship , where τ is the relaxation time, E is the activation\na\nenergy, k is the Boltzmann constant, and T is the temperature. As the temperature increases, the\n𝜏 = 𝜏0𝑒𝑥𝑝(𝐸𝑎⁄𝑘𝑇)\ndipoles acquire sufficient thermal energy to overcome the potential barriers, resulting in faster\npolarization dynamics (smaller τ). Consequently, the dielectric peak moves to higher\ntemperature with increasing frequency in typical relaxor systems.\nHowever, in the CuInP (S Se ) crystals, the substitution Se to S introduces defect states\n2 1−x x 6\nand additional carriers, leading to enhanced electrical conductivity. At low frequencies of 300\nHz, the carriers have sufficient time to respond to the external field, and space-charge\naccumulation at grain or domain boundaries gives rise to pronounced Maxwell-Wagner-type\ninterfacial polarization. In this case, the overall dielectric response arises from the superposition\nof dipolar and space-charge polarizations, and the relaxation process corresponding to the\nmaximum dielectric constant occurs at relatively higher temperatures. As the frequency\nincreases up to 3 MHz, the space-charge contribution is gradually suppressed since the charge\ncarriers can no longer follow the rapid electric field. The dielectric response then becomes\ndominated by intrinsic dipolar polarization, whose relaxation time corresponds to lower\ntemperature. Therefore, the dielectric peak shifts toward a lower temperature with\nincreasing frequency.\nThis anomalous shift indicates that the Se-doped CuInP (S Se ) is not a typical\n2 1−x x 6\n7\nrelaxor ferroelectric, but rather a composite relaxation system arising from the\ncoexistence of conductive space-charge and dipolar polarization. And the similar\nphenomena have been reported in Fig. R6 (Science 385, 204–209 (2024)) and Fig. R7 (Nat.\nMater. 19, 1169–1174 (2020)).\n[Figure Redacted]\nFig. R6 (Adopted from Science 385, 204–209 (2024)) Temperature-dependent dielectric\npermittivity and loss tangent of the (A) PNR-type, (B) PS-type, and (C) IPS-type BMT-ST-\nbased RFE films and their dielectric permittivity variation around T .\nm\n[Figure Redacted]\nFig. R7 (Adopted from Nat. Mater. 19, 1169–1174 (2020)) Extend Data Fig. 7a, The red line\nis a fit of the dielectric constant (solid black symbols) by the Vogel-Fulcher law. Extend Data\nFig. 8, Conductivity at different temperatures and frequencies. The gray region indicates the\ndominant region of ionic conductivity.\nTo address comment, we have supplied this explanation into the supplementary Note 3\nas follows: “In the Se-doped CuInP (S Se ) crystals, the increased electrical conductivity\n2 1–x x 6\ncaused by Se substitution introduces significant Maxwell–Wagner polarization at low\nfrequencies. At the low frequency, the charge carriers can accumulate at internal interfaces or\ndefects, producing an additional space-charge polarizations. In this case, the overall dielectric\nresponse arises from the superposition of dipolar and space-charge polarizations, and the\nrelaxation process corresponding to the maximum dielectric constant occurs at relatively\nhigher temperatures. As frequency increases, this space-charge contribution is gradually\nsuppressed because carriers cannot follow the rapid electric field. The dielectric response then\nbecomes dominated by faster dipolar relaxation, which maximum dielectric constant occurs at\nlower temperatures. Thus, the dielectric peak shifts toward lower temperature with increasing\nfrequency.”\niii. Is it possible to estimate the phase fraction of the trigonal phase with increase in the Se\ndoping concentration?\nResponse:\nWe sincerely thank you for your insightful suggestion. The variation of the trigonal and\nmonoclinic phases was mainly analyzed by XRD and Raman. As shown in the XRD spectra\n(Fig. R8), the (002) diffraction peak gradually shifts to higher angles, indicating a reduced\ninterlayer spacing along the c-axis, which is associated with the emergence of the trigonal phase\nwith smaller layer distances. Meanwhile, the intensities of the (004), (006), and (008) peaks\ngradually decrease, suggesting that the uniformly distributed Se promotes a higher proportion\n8\nof trigonal phase and disrupts the long-range order of crystal, which is consistent with the polar\nnanoregions. In the Raman spectra (Fig. R9), the characteristic peak at 315 cm⁻1 (corresponding\nto the Cu hopping mode) progressively weakens, while the peak near 346 cm⁻1 (assigned to the\nvibration mode of the trigonal phase) is enhanced with higher Se content. These results\nconsistently indicate that the fraction of the trigonal phase increases with Se doping.\nHowever, these results can only provide the relative phase evolution rather than precise\nvolume fractions. This is because the Raman peak intensity is not linearly correlated with the\nphase fraction, and the XRD intensity is also affected by crystallite size, amorphous regions,\nand dislocation density. Furthermore, although STEM imaging can directly reveal the local\ndistribution of the two phases, its limited field of view cannot accurately represent the actual\nphase composition within the crystal. Therefore, combined these factors, the precise phase\nfractions cannot be obtained, while an increasing trend of the trigonal phase fraction could be\nverified by the XRD and Raman results.\nFig. R8 (Updated in Fig. 2a in Manuscript) XRD patterns of Se content increasing from 0%\nto 6.9%.\n2 1-x x 6\nwith Se content increasing from 0% to 6.9%.\niv. The authors highlight the importance of the strain gradients under Se substitution towards\nthe formation of the trigonal phase. My question in this regard is related to point 1 above. How\nfeasible will it be to generate such strains during the growth of a just a few monolayer thick\nfilms due to the differences in the mechanical boundary conditions between crystal growth vs\nfilm growth?\nResponse:\nThanks for your professional question. The 2D crystals were primarily synthesized via the\nchemical vapor transport (CVT) method, which enables the formation of high-quality single\ncrystals under thermodynamic equilibrium growth conditions. In this case, the strain was\nintroduced mainly through compositional tuning via doping. Specifically, in the\nCuInP (S Se ) system, the partial substitution of S by Se atoms, which have larger ionic radius\n2 1−x x\nand lower electronegativity, causes the increased lattice instability, thereby introducing a strain\ngradient within the crystal. This is fundamentally different from thin-film growth methods,\nwhere strain is often modulated through the substrates, deposition conditions, or annealing\nprocesses. Furthermore, since the substitution occurs within the layers, the resulting strain\ngradient remains preserved even after the crystal is exfoliated into few-layer flakes. This also\nresponse the point 1 of your review, by explaining the persistence of relaxor ferroelectric in the\nfew-layers crystals.\nv. What is the physical and technological significance of the relaxor diffuseness factor with Se\ndoping of x=1.4% being 1.45 mean? At x=4.2%, and the diffuseness factor of 1.9 seems to\n10\nsuggest a much stronger relaxor-like behavior, but this phase is stabilized only at low\ntemperatures but for most practical technological applications, the temperature is usually RT or\nabove.\nResponse:\nWe really appreciate the perceptive question. The relaxor diffuseness factor (γ), as derived\nfrom the modified Curie-Weiss law, typically ranges between 1 and 2, where γ = 1 corresponds\nto a typical ferroelectric behavior, and γ = 2 represents an ideal relaxor ferroelectric. In\nthis study, the crystal with 2.6% (corresponding to 1.4% Se content in previous\nmanuscript) Se content exhibits a γ value of 1.45, indicating the formation of a relaxor\nferroelectric phase with a diffuse phase transition, where local polar nanoregions are\nappeared. As the Se content increases to 6.9% (corresponding to 4.2% Se content in previous\nmanuscript), γ further increases to 1.9, reflecting a more pronounced relaxor behavior\ncharacterized by typical features of a diffuse phase transition.\nThe relaxor diffuseness factor γ is derived from the modified Curie-Weiss law:\n1 1\n𝜀𝑟 −𝜀𝑟,𝑚 =\n, in which ɛ is the permittivity at the temperature T, T is the temperature at which ɛ\nr m r\n𝛾\n(𝑇−𝑇𝑚)\nreaches the maximum ɛ , and C is a constant. In a typical ferroelectric, the low-temperature\n𝐶 r,m\nferroelectric phase transforms into the high-temperature paraelectric phase at the Curie\ntemperature (T ), where the dielectric permittivity reaches maximum. Above this temperature,\nc\nthe spontaneous polarization vanishes rapidly, leading to a sharp drop in permittivity. Such\ntransition follows the standard Curie-Weiss law, corresponding to a γ ≈ 1.\nIn contrast, the relaxor ferroelectric exhibits a typical diffuse phase transition, where the\nCuire temperature is not fixed value but distributed over a range, forming the Curie temperature\nregion. This arises from the presence of polar nanoregions and polarization diminishes\ngradually with temperature over the T . Consequently, the dielectric peak becomes broader and\nc\nthe Tc exhibits a frequency dependent shift. This phase transition can be described by the\nmodified Curie-Weiss law.\nTherefore, the relaxor diffuseness factor γ reflects the breadth of the distribution of phase\ntransition temperatures. When γ approaches 1, the system exhibits a typical ferroelectric,\nobeying the standard Curie-Weiss law with fixed Cuire temperature and sharp dielectric peak\nthat is frequency independent. As γ increase toward 2, the system gradually transforms to\nrelaxor ferroelectric, which is characterized by diffuse phase transition, the dielectric peak\nbroadening and frequency dependence arising from the emergence of polar nanoregions.\nDue to the presence of diffuse phase transition in relaxor ferroelectric, the dielectric peak\narounds the Curie temperature becomes broadened, and the polarization decreases gradually\n11\nabove Tc. This behavior enables a polarization response over a wide temperature range. Even\nthough the Tc shifts down to room temperature with increasing Se content, the partial\npolarization could still be remained above the room temperature in the relaxor ferroelectric\ncrystal. For instance, the crystal with Se content of 6.9% exhibits a stable dielectric constant of\n60 over the wide temperature range of –100 °C to 50 °C, whereas the undoped crystal (x = 0%)\nwith dielectric constant of only 14 within a narrow range from –100 °C to 20 °C (Fig. R10).\nThis provide clearly evidences that relaxor ferroelectric exhibits a broader temperature window\nfor polarization response or potential applications.\nFig. R10 (In Supplementary Fig. 4b in Supporting Information) Temperature dependent of\ndielectric permittivity and tan (δ) with Se content increase from 0% to 6.9%.\nTo address this comment, we add “where the ε is the maximum dielectric constant at T ,\nm m\nand C is the Curie-like constant, where the γ = 1 corresponds to a typical ferroelectric behavior\nand γ = 2 represents an ideal relaxor ferroelectric. The γ values increase progressively from\n1.11 for x = 0% to 1.45 for x = 2.6%, and 1.9 for x = 6.9%, demonstrating the enhanced relaxor\nproperty (Fig. 3f). And the γ of 1.45 at x = 2.6% indicates the initial transition from FE toward\nto RFE, which corresponding to the emergence of phase coexistence structure.” on page 6 of\nthe revised manuscript. And “This enhancement is attributed to the pronounced relaxor\nbehavior, which enables the polarization to gradually vanish even above the T . And the\nm\nincrease of dielectric permittivity is associated with the local polar regions, which allows larger\npolarization response at lower electric field.” on page 7 of the revised manuscript.\nvi. Can the authors elaborate more on why they considered relaxors for potential memristor\napplications? Relaxors are known to have lower remanent polarization which would mean that\nthey are more susceptible to retention loss which would seem detrimental for any type of\nmemory applications. It might have been more useful to quantify the electrostrictive\ncoefficients or measure the strain in such relaxor systems.\nResponse:\nThank you for the opportunity to clarify this point. As you point out, the relaxor have low\nretention time, which means it is not suitable for long-term information storage. However,\nbenefitted from their excellent multi-state memory feature arise from gradual switching process\nof polar nanoregions, we note that they are quite compatible to the recently-developed Artificial\nNeutral Network (ANN) hardware. For such a neuromorphic computing, it is constructed by\nthe large array of numerous crossbar memory cells, which is similar to this two-terminal\nmemristor devices. Most importantly, such an application scenario requires each cell to be able\nto fast adjust their conductance between multiple levels, in order to endow the “learning” ability\nto this hardware system. Therefore, the long-term memory ability is not so necessary here,\ninstead, the multi-level feature of our devices is more emphasized in such scenario.\nFollowing your suggestion, the piezoelectric strain was measured at the resonance\nfrequency using PFM, and the results show that the piezoelectric strain decreases as the phase\ntransforms from the ferroelectric to relaxor ferroelectric (Fig. R11). This decrease is associated\nwith the presence of a lower-symmetry trigonal phase after Se doping, which reduces the\npolarization.\nFig. R11 Piezoelectric strain of crystal with Se content of 0%, 2.6% and 6.9% under drive\nvoltage of 1V.\nTo address this comment, we add “Therefore, the formation of nanodomains in relaxor\nferroelectric enables multiple resistance state modulation under lower operating voltages,\nmaking it well-suited for the memristor devices used in neuromorphic computing.” on page 7\nof revised manuscript.\n13\n3. Other comments:\ni. Please use a consistent system of labeling the supplementary figures either calling them\nSupplementary Figure x or Figure Sx.\nii. Fig. 3g: the figure says the transistion is from AFE to FE but the label says it shows transition\nfrom AFE to RFE?\niii. Fig.2c: labels of atoms are missing\nResponse:\nTo address this comment:\ni. We have standardized the labeling of all supplementary figures throughout the manuscript to\n“Supplementary Fig. x” and have carefully checked to ensure consistency.\nii. Thank you for pointing out the inconsistency in Fig. 3g. We have corrected the figure label\nso that they consistently indicate the transition from AFE to RFE in Fig. R12.\niii. We have added the missing atomic labels in Fig. R13 to make the structural assignment clear.\nFig. R12 (Updated in Fig.3g in Manuscript) Phase transition from antiferroelectric to relaxor\nferroelectric at freezing temperature (T) with Se content of 2.6%.\nf\ncrystal at Se content of 2.6% showing the monoclinic (M), trigonal (T) phase and edge\ndislocation, scale bar: 5 nm. The edge dislocation region is marked with orange dashed lines.\nThe diffraction patterns show monoclinic and trigonal phases. Geometric phase analysis (GPA)\nwas applied to extract the strain distribution at the edge dislocation regions.\nAll corrected figures and captions are included in the revised manuscript. We appreciate\nthe reviewer’s careful reviewing, which helped us improve the clarity and accuracy of the\npresentation.\nThe manuscript by Yang et al. is entitled “Phase engineering of relaxor ferroelectricity in van\nder Waals crystal” (id: NCOMMS-25-70669-T), which realizes the relaxor ferroelectricity in\nlow-dimensional materials by controlling phases. The phase engineering is realized by partial\nsubstitution of S with Se in CuInP (S Se ) and establishes competing coexistence of\n2 1-x x 6\nmonoclinic and trigonal phases. This approach reveals a sequential transformation of FE to RFE\nand SPE phases at room temperature in the 2D materials. Through experimental techniques like\nXRD, SAED, TEM images, SHG signal, and PFM, they tried to establish their goal. The data\nquality seems sufficient to support the main conclusions. They use DFT to compare the energy\nbarrier of phase transition and polarization of the two pristine monoclinic and trigonal phases.\nMost of the manuscript is well written; however, it would benefit from a careful revision to\neliminate repetitions, correct typos, and fix mislabeled references. Before recommending\npublication in Nature Communications, the authors should consider the following points.\nResponse:\nWe sincerely appreciate your careful review and thoughtful evaluation of our manuscript.\nMeanwhile, we are grateful for your encouraging comments on our work as “Most of the\nmanuscript is well written;”. In response to these valuable comments, we have meticulously\nrevised the manuscript to improve the quality of the study. We believe that these revisions have\nmade our work clearer, more comprehensive, and scientifically more robust.\ni) The authors claim that they have realized the relaxor ferroelectricity in 2D materials for the\nfirst time. However, there is a report that previously indicated the realization of relaxor\nferroelectricity in Se-doped 2D PbTe by Gong et al. (Strain-induced two-dimensional relaxor\nferroelectrics in Se-doped PbTe; DOI: https://doi.org/10.1103/PhysRevB.109.054117). It is not\nclear why the authors insist on this statement. Either the authors should be clearer with their\nstatement or remove such a statement.\nResponse:\nWe really thank you for raising this important point. The previously study, “Strain-induced\ntwo-dimensional relaxor ferroelectrics in Se-doped PbTe”, offers only theoretical predictions,\nwhile our work provides the first experiments evidence. The reported relaxor ferroelectric\nin PbTe was predicted computationally but has not yet been experimentally validated. In\ncontrast, our work presents the first experimental realization of relaxor ferroelectricity in a two-\ndimensional material, providing a solid foundation for further exploration of 2D relaxor\nferroelectric properties. To avoid this misunderstanding, we have specifically emphasized that\nthis is the first time relaxor ferroelectricity has been experimentally achieved in a synthesized\n2D crystal.\nin low-dimensional materials for the first time.” on page 1 of revised manuscript.\nii) As mentioned, the XRD diffraction peaks at 13.4° correspond to the (002) indicates the 2D\ncrystal stacks along the c-axis. According to the Authors, broadening the (002) diffraction peak\nsuggests that atomic substitution induces grain size refinement, interlayer dislocations, and\namorphous phases for Se content of x = 1.4%. However, when the concentration of Se further\nincreased to x = 2.5% and x = 4.2%, the (002) peak again became sharp in nature. The reason\nfor this anomaly is not discussed in the manuscript. Authors are advised to discuss this to clear\nthe anomaly.\nResponse:\nWe greatly appreciate the careful review and valuable suggestion. The sharpening of the\ndiffraction peaks at Se contents above 2.6% (corresponding to 1.4% Se content in previous\nmanuscript) is mainly attributed to increased lattice ordering. As the Se content increases\nto 2.6%, the local substitution of Se produces a gradient strain large enough to distort the lattice,\nwhich disrupts the long-range stacking order and leads to edge dislocation as well as the\ncoexistence of trigonal and monoclinic phases. This structure heterogeneity results in the\nbroadening of the (002) diffraction peak. When the Se content is further increased, the\nsubstituted Se atoms become more uniformly distributed, effectively releasing the accumulated\nstrain, reducing dislocation density, and improving the lattice order. Consequently, a sharper\n(002) diffraction peak appears at x = 4.2% and 6.9% (corresponding to 2.5% and 4.2% Se\ncontent in previous manuscript). Therefore, x = 2.6% represents a critical threshold in Se site\noccupancy, marking the transition in dislocation within the crystal during the Se substitution.\nTo address this comment, we add “As the Se content surpass 2.6%, the more uniform\ndistribution of Se effectively releasing the accumulated strain and improving the local stacking\norder, resulting in the recovery of sharper (002) diffraction peak.” on page 4 of revised\nmanuscript.\niii) Since these are the multilayers of CuInP2(S Se ) . There is a high possibility of having\n1-x x 6\ndifferent concentrations of Se substitution in the different layers, leading to differences in the\nlattice parameter. Also, as the coupling between the layers is expected to be weak, a moiré\npattern should emerge; did the authors check for it with STM? If they have observed that, then\nwhat is the role of that varying interlayer coupling for this phase transformation?\nThank you for your insightful suggestion. To reveal the moiré pattern in atomic-resolution\nimages, the HAADF-STEM was performed along (001) crystal direction, and moiré pattern\nwas not observed. The explanation is that the lattice parameters of the trigonal and monoclinic\nphases is too similar to produce the moiré pattern.\nBased on the simulated diffraction patterns and atomic lattices of monoclinic and trigonal\nphase (Fig. R14)3, the typical monoclinic phase was identified in Fig. R15a. The additional\ndiffraction spots observed in Fig. R15b indicate the nucleation of a trigonal phase, while the\nbrighter trigonal reflections in Fig. R15c suggest an increased fraction of the trigonal phase.\nThe gradual transition from the monoclinic to the coexisted-phase is also evident in the real-\nspace atomic images. And the moiré pattern was not observed in the STEM images, and the\nmoiré-related satellite reflections were not measured in the corresponding FFT patterns,\nindicating the absence of moiré pattern in the coexisted-phase region.\nThis is attributed to the fact that the monoclinic and trigonal phases possess nearly identical\nin-plane lattice parameters. In the monoclinic phase, the lattice constants are a = 6.09 Å and b\n= 10.56 Å, with an in-plane angle of 90°. In contrast, the trigonal phase has a = b = 6.09 Å with\nan in-plane angle of 120°. The similarity of the in-plane lattice parameters reaches 99.8%,\nwhich is insufficient to induce fixed lattice mismatch. Therefore, the moiré patterns are not\nobserved at the coexisted phases.\nFig. R14 Simulated monoclinic and trigonal structures of CuInP S . a,b The structural of the\n2 6\n18\nmonoclinic and trigonal phases viewed along the c-axis ((001) direction). c,d Simulated\nHAADF-STEM atomic images of the monoclinic and trigonal phases. e,f Simulated electron\ndiffraction patterns for the monoclinic and trigonal phases.\nFig. R15 (Added in Supplementary Fig. 5 in Supporting Information) HAADF-STEM\nimages and corresponding fast Fourier transform (FFT) patterns of CuInP (S Se ) crystals\n2 1−x x 6\nacquired along the (001) direction. a, Monoclinic phase showing a well-ordered lattice and FFT\ndiffraction pattern. Scale bar, 1 nm. b, Mixed-phase region with the coexistence of monoclinic\nand trigonal phases, as indicated by additional trigonal diffraction spots in the FFT. Scale bar,\n1 nm. c, Mixed-phase region with a higher fraction of the trigonal phase, evidenced by brighter\ntrigonal reflections in the FFT and further lattice distortion in the real-space image. Scale bar,\n1 nm.\nReferences:\n3. Deng, J. et al. Thickness-dependent in-plane polarization and structural phase transition in\nvan der Waals ferroelectric CuInP S . Small 16, 1904529 (2020).\n2 6\nTo address this comment, we have updated the Supplementary Fig. 5 by including the\nHADDF-STEM images of the crystal along the (001) direction, to better illustrate the evolution\nof atomic-scale contrast with the increased trigonal phase fraction.\niv) Using DFT, the Authors tried to support the easy transition along FE-RFE-SPE. However,\nthe results of Fig. 4(c) and 4(d) are for pristine CuInP Se and CuInP S . Furthermore, this\n2 6 2 6\nresult is derived by connecting the P_down to the P_up structure. The easy transition along FE-\nRFE-SPE, using these pristine calculations, does not make strong evidence for the experimental\nsupport. The calculations using CuInP2(S Se ) structure will give more realistic evidence.\n1-x x 6\nResponse:\nWe appreciate for the insightful suggestion. The ferroelectric behavior of the mixed-phase\ncrystal is primarily governed by the coexistence of the monoclinic (CuInP S ) and trigonal\n2 6\nphases (CuInP Se ), where the appearance of the coexisted phase induces locally polarized\n2 6\nregions. Therefore, we calculated the polarization characteristics of both the monoclinic and\ntrigonal phases using DFT. The results reveal that with increasing trigonal phase fraction, the\nsystem gradually evolves from FE to RFE and SPE, accompanied by a reduction in both\npolarization strength and coercive field. This trend is consistent with the lower energy barrier\nand smaller polarization obtained for the trigonal phase in our calculations.\nWe fully agree that DFT calculations based on the CuInP (S₁₋ Se ) structure would\n2 x x 6\nprovide a more realistic representation of the local atomic environment. After discussed with\nexperts on DFT, we found such modeling requires extremely large atomic models and\ncomputing time, which is extremely challenging and computationally not practical. Since\nthe Se concentration in the mixed crystal is relatively low, constructing a realistic random\nsubstitution model requires very large supercells to capture the stochastic distribution of Se\natoms. Each disordered configuration must be fully relaxed to account for the local strain and\ninterlayer dislocations, which dramatically increases the computational cost. Moreover,\nexploring the FE-RFE-SPE transition pathway would require performing polarization reversal\ncalculations for each configuration, rendering the total DFT workload computationally\nprohibitive.\nTherefore, direct modeling of the CuInP (S₁₋ Se ) through first-principles calculations is\n2 x x 6\nnot feasible at the current stage. Instead, we analyzed the polarization switching behavior of the\nindividual monoclinic and trigonal phases to elucidate the FE-RFE-SPE transition in the mixed-\nphase crystal. This theoretical trend is consistent with the experimentally observed reduction in\ncoercive field and polarization strength upon Se incorporation.\nv) The Authors mentioned that there is a decrease in switching voltages for domain reversal\nwith increasing Se content (Fig. 4b). This claim was valid and consistent from x = 0% and x =\n2.5%. However, for x = 4.2%, the switching voltage is higher than x = 1.4% and 2.5% on the\npositive side. Please clear this anomalous observation.\nResponse:\nWe sincerely appreciate your professional question. During the phase loop measurement,\n20\nthe switching voltage is influenced by both the spontaneous polarization and contact potential\ndifferences between the conductive tip and the sample surface. These factors can induce an\ninternal built-in electric field, leading to asymmetric energy barriers for polarization reversal.\nAs a result, the measured phase loops exhibit an asymmetric shift with respect to zero bias, as\nobserved in Fig. 4b. However, this does not affect the conclusion of the switch voltage\ndecrease, as the effective coercive voltage (Vc) can be defined as , where Vs+\n+ −\n|𝑉𝑠 |+|𝑉𝑠 |\nand Vs– denote the positive and negative switching voltages, respectively. Based on this\n𝑉𝑐 = 2\ndefinition, the coercive voltage is found to gradually decrease from 7.14 V to 2.8 V as the Se\ncontent increases from x = 0% to x = 6.9% (corresponding to 4.2% Se content in previous\nmanuscript).\nMinor comment:\ni) The texts in SI from line numbers 66 to 106 have an issue. Also, after line number 106, line\nnumber 265 has started.\nResponse:\nTo address this comment, we carefully inspected the Supporting Information and\nidentified a formatting error that caused disruption in the text block and line numbering between\nlines 66-106 and line 265. The issue has been corrected, and the line-numbering continuity has\nbeen restored. The revised Supporting Information has been reuploaded accordingly.\nIn this manuscript, the authors addressed an interesting work on phase engineering of relaxor\nferroelectricity in van der Waals crystals. Authors demonstrated a mixed-phase 2D CuInP (S\n2 1-\nSe ) single crystals where tuning the ratio of monoclinic and trigonal phases induced\nx x 6\nsuperlattice structures and heterogeneous local polar states, enabling transitions from\nferroelectric to relaxor, superparaelectric, or paraelectric behavior through a flattened\nfree-energy profile. Structural, compositional, electrical, and atomic scale measurements were\nperformed using various characterization techniques. Temperature-dependent dielectric and\nferroelectric measurements confirmed relaxor behavior via hysteresis-free polarization and\nfrequency dispersion of permittivity, while PFM revealed domain size and coercive voltage\nevolution. Integration into memristor devices further showed enhanced analog resistance states\nand low operating voltages, highlighting potential for electronic applications. Although the\nmanuscript was interesting. However, I regret that the explanation and mechanism of the\nmanuscript were not convincing to support many of their major claims. The mechanism of the\nstudy remained confusing due to incomplete quantitative evidence and limited analysis. In its\ncurrent form, the manuscript cannot be considered for publication until all the raised concerns\nwere properly addressed.\nResponse:\nWe would like to express our true appreciate for your time and effort in reviewing our\nmanuscript and thorough evaluation. We are truly grateful for pointing out the important\nprofessional comments. We carefully examined the questions, and conducted additional\nexperiments to strength the evidence for the mechanism understanding. Quantitative evidence\nand detailed analyses have been added to the revised manuscript. We hope that the revised\nversion strengths the mechanism analysis and provides sufficient discussion to meet the high\nstandards of this journal.\nMajor comments\n1) A central concern was the issue of composition control and quantification. The reason for\nselecting a non-uniform substitution trend of the Se ratio into the host lattice was unclear.\nTypically, to establish compositional trends and phase boundaries, systematic increment doping\nstrategies should be employed. The current approach raised questions about how the authors\ncontrolled and verified the exact Se content, particularly at the extremely low substitution\n(around 1.4%) levels reported. Without clear justification and compositional control, the\nscientific basis of the correlation between composition and physical behavior was unclear.\nResponse:\nWe sincerely appreciate for the insightful question. To capture the composition\n22\ndependent phase transition, the initial Se ratios in the precursor mixtures were set as 0%,\n1%, 3%, 5%, and 8%. Unlike in ceramic systems, the elemental composition is identical to\nthe elemental ratios in the precursors, in this CVT grown van der Waals system, a small portion\nof Se elements could be excluded from the lattice, resulting in slight reduction of actual Se\nincorporation. Through experiment measurements, we concluded that the actual composition\nratio after growth is approximately 0.2% lower than the value set in the precursor. This is\nattributed to the differences in vapor-transport efficiency and minor elemental weighing errors\nduring the growth process. Therefore, the actual Se content in the grown crystals does not\nexactly match the nominal values.\nTo obtain the accurate Se content, the EPMA measurements was performed. The\nrefined Se contents for the nominal 1%, 3%, 5%, and 8% samples were determined to be\napproximately 0.8%, 2.6%, 4.2%, and 6.9%, respectively (Fig. R16). The actual Se\nincorporation remains consistent trend with the initial precursor ratios.\nWhile the compositional increase was not perfectly linear, the selected Se concentrations\neffectively captured the critical transition in ferroelectric behavior, reflecting the coupled\nevolution between structural modulation and property. The systematic evolution of XRD\ndiffraction peaks, Raman vibration modes, and SHG responses with increasing Se content\nprovides strong indirect evidence that the Se incorporation was well controlled and gradually\nincreased. These results together provide a reliable basis for correlating structure and property\nrelations in the crystal.\nFig. R16 Comparison between the actual Se content and the initial elemental ratios in the\nprecursors for CuInP2(S Se ) crystals.\n1-x x 6\nTo address this comment, We add “High pure precursor materials, including Cu (99.99%,\nSigma-Aldrich), In (99.99%, Sigma-Aldrich), P (99.999%, Sigma-Aldrich), S (99.98%, Alfa\nAesar) and Se (99.99%, Sigma-Aldrich) were mixed in quartz tube with a nominal\n23\nstoichiometric proportion of 1:1:2:6:0, 1:1:2:5.982:0.018, 1:1:2:5.916:0.084, 1:1:2:5.85:0.15,\n1:1:2:5.748:0.252, corresponding to the Se content 0%, 1%, 3%, 5% and 8%, respectively. This\nslightly deviation is ascribed to the variations in vapor transport efficiency and minor weighing\ninaccuracies of the elemental precursors during the CVT growth process.” on page 14 of revised\nmanuscript.\n2) The authors synthesized a series of CuInP (S Se ) crystals with the increasing Se content\n2 1-x x 6\nof 0%, 0.3%, 1.4%, 2.5%, and 4.2%. Based on the formula, the Se should be homogenously\nsubstituted rather than the formation of a composite single crystal, where a mixed region of\nmonoclinic and trigonal would exist with shared boundaries, where two regions apply\ninterfacial strain. In composites, lattice mismatch at interfaces could generate strain. The\nauthors should discuss why uniform substitution was pursued over this alternative and whether\nthe substitution was homogeneous or whether the authors synthesized the composite single\ncrystals.\nResponse:\nWe highly appreciate the professional question. For the CuInP (S Se ) crystals, the\n2 1-x x 6\nsubstituted Se atoms tend to randomly occupy the S sites leading to macroscopically\nhomogeneously distribution, while exhibiting local inhomogeneous substitution at the\nmicroscopic level. This macroscopically uniform distribution of Se is further verified by the\nEPMA elemental mapping (Fig. R17). The local inhomogeneous substitution is closely related\nto the emergence of dislocations. When the Se concentrations below 2.6% (corresponding to\n1.4% Se content in previous manuscript), the substituted Se atoms are randomly and sparsely\ndistributed within the 2D layers, and the induced point defects are insufficient to distort the\nlattice, maintaining a homogeneous single-phase structure. As the Se content increases to 2.6%,\nlocal Se enrichment near defect sites induces sufficient strain to distort the lattice, leading to\nthe formation of dislocation and trigonal phase. Upon further increasing the Se content beyond\n2.6%, the uniformly distributed Se atoms effectively release the accumulated strain and\ndecreases dislocation density, result in local heterogeneity and the coexistence of monoclinic\nand trigonal phases. Therefore, this uniform substitution was favored over this alternative,\nand the composite single crystal did not form because of the low substitution. This type of\nnon-stoichiometric substitution-induced local phase coexisted is commonly observed as a\nstructural modulation strategy in oxide systems, as reported in previous studies 4-7.\nwith a Se content of 2.6%. Scale bar, 50 µm. b, EPMA elemental mapping of Se showing its\nspatial distribution across the crystal surface. Scale bar, 50 µm.\nReferences:\n4. Pan, H. et al. Ultrahigh–energy density lead-free dielectric films via polymorphic\nnanodomain design. Science 365, 578-582 (2019).\n5. Chen, L. et al. Giant energy-storage density with ultrahigh efficiency in lead-free relaxors\nvia high-entropy design. Nat. Commun. 13, 3089 (2022).\n6. Shu L, et al. Partitioning polar-slush strategy in relaxors leads to large energy-storage\ncapability. Science 385, 204-209 (2024).\n7. Du F, et al. Giant electrocaloric effect in high-polar-entropy perovskite oxides. Nature, 640,\n924-930 (2025).\nTo address this comment, we add “During the CVT growth, the Se atoms incorporate\ninto the lattice under near-equilibrium conditions, leading to macroscopically homogeneously\ndistribution at S sites. The monoclinic and trigonal phases can be determined by the angles\nbetween a- and c-axis and the interlayer distance along the (010) plane (Supplementary Note\n4, and Supplementary Fig. 7). For Se content below 2.6%, the substituting Se atoms are\nrandomly and sparsely distributed within the 2D layers, and the induced point defects are\ninsufficient to distort the lattice, maintaining a homogeneous single-phase structure. And the\nindependent monoclinic phase was observed, characterized by the angle of 98.9° between the\na and c axis and clearly vdW gap (Supplementary Fig.5 and Fig. 8). With the Se content increase\nto 2.6%, local Se enrichment near defect sites induces sufficient strain to distort the lattice,\nleading to the formation of dislocation and trigonal phase.” on page 4 of revised manuscript.\n3) In Supplementary Table 1, the atomic percentage values of Cu, S, and Se were inconsistent,\nwhich challenges the compositional analysis. EDS alone was insufficient for quantitative\ndetermination, especially at such low substitution levels. The authors should validate their\nresults with more precise techniques, such as the electron probe microanalyzer (EPMA)\ntechnique [Phys. Rev. Materials 3, 084405 (2019)].\nResponse：\nThank you for your valuable suggestion for the composition analysis. Following the\ncomment, we have conducted the EPMA measurement to accurately quantify the Se\nsubstitution level. For the nominal Se contents of 1%, 3%, 5%, and 8%, the accurately\nquantified Se concentrations were determined to be approximately 0.8%, 2.6%, 4.2%, and 6.9%,\nrespectively (Table R1). And each sample of the Se content was determined by averaging four\nmeasurement points to ensure statistical reliability. These results have already been discussed\nin the responses to your comments: Qestion1 and Qestion2.\nTable R1. (Updated in Supplementary Table1 in Supporting Information) Calculated\natomic ratios obtained from EPMA measurements for crystal with gradually increased Se\ncontent.\nTo address this comment, we have updated Supplementary Table1 include the accurate\nSe contents of 0%, 0.8%, 2.6%, 4.2%, and 6.9% in the revised manuscript, as obtained from\nEPMA measurements.\n4) Phase identification was another weak point of the manuscript. The authors showed the\ncoexistence of monoclinic and trigonal phases, and the proposed superlattice formation,\nqualitatively. But there was no quantitative presentation, like reporting of lattice parameters or\nphase fractions.\nResponse:\nIn the monoclinic and trigonal phases, the differences in lattice parameters give rise to\ncorresponding variations in the crystal structure. As previously reported3,8, the monoclinic\nphase exhibits lattice constants of a = 6.09 Å, b = 10.56 Å, and c = 13.17 Å, whereas the\ntrigonal phase shows a = b = 6.09 Å and c = 13.00 Å. Another notable distinction is that the\nangle between the a and c axes is approximately 98.9° in the monoclinic phase, while it is\n90° in the trigonal phase (Fig. R18). These crystallographic differences allow us to distinguish\nthe two phases in atomic-resolution STEM images obtained from different samples by\nanalyzing the lattice orientation within the ac plane (Fig. R19).\nthe trigonal structure, the Cu atoms shift slightly inner layer and the S atoms move outward,\nleading to a reduction in the van der Waals interlayer spacing (Fig. R19). As a result, the\ninterlayer distance in the trigonal phase (3.116 Å) is smaller than that in the monoclinic phase\n(3.244 Å). Consistently, the c-axis lattice constant of the monoclinic phase is larger than that of\nthe trigonal phase. This trend is corroborated by the XRD patterns, in which the (002)\ndiffraction peak progressively shifts to higher angles with increasing Se concentration,\nconsistent with a reduction of the average interlayer spacing from 13.2 Å to 13.1 Å due to the\nenhanced proportion of the trigonal phase (Fig. R20).\nFig. R18 (Updated in Supplementary Fig. 7 in Supporting Information) Crystal structure\nof CuInP S with monoclinic and trigonal phase. a, Monoclinic phase structure along the c, a\n2 6\nand b axes, respectively. The angle between a and c is 98.9°. b, Trigonal phase structure along\nthe c, a and b axes, respectively. The angle between a and c is 90°.\nFig. R19 The interlayer distance of monoclinic and trigonal phase. The calculated d is 3.116\n1\nÅ and the d is 3.244 Å.\n2\ncontent increasing from 0% to 6.9%.\nReferences:\n8. Li, T. et al. Realization of sextuple polarization states and interstate switching in\nantiferroelectric CuInP S . Nat. Commun. 15, 2653 (2024).\n2 6\nTo address this comment, we have updated Supplementary Fig. 7 of the crystal structure\nwith monoclinic and trigonal phase. The corresponding lattice parameters was supplied as “In\naddition, the structure difference between the two phases enables the phase identification. The\nmonoclinic phase exhibits lattice constants of a = 6.09 Å, b = 10.56 Å, and c = 13.17 Å, whereas\nthe trigonal phase shows a = b = 6.09 Å and c = 13.00 Å10,11. Another notable distinction is\nthat the angle between the a and c axes is approximately 98.9° in the monoclinic phase, while\nit is 90° in the trigonal phase (Supplementary Fig. 7). These crystallographic differences allow\nus to distinguish the two phases in atomic-resolution STEM images by analyzing the lattice\norientation within the ac plane. In the HAADF-STEM images along the (010) plane, the\nuniformly distributed monoclinic phase is clearly observed at Se content of 0% (Supplementary\nFig. 8).”\n5) Based on the STEM results, the authors attributed the formation of edge dislocations to the\nSe substitution, but the doping levels were too low to justify such noticeable defect generation.\ngrowth-related stresses, must be considered.\nResponse:\nThank you very much for providing these possible explanations of the origin of edge\ndislocations. We fully agree with your view that strain gradients at grain boundaries are\ncritical factor contributing to detect formation. Given the larger ionic radius and weaker\nelectronegativity of Se atoms compared to S, Se substitution preferentially occurs at grain\nboundaries. This leads to the formation of local Se-rich regions, which introduce noticeable\nlattice distortions and generate substantial local strain. To minimize the elastic energy, the\ndislocations are formed in these regions. This mechanism is supported by the STEM-GPA strain\nmapping results, which clearly show pronounced shear stress around the dislocation cores (Fig.\nR21). These localized strain fields further induce interlayer sliding within the layered crystal\nstructure, leading to the coexistence of monoclinic and trigonal phases. Therefore, the small\nSe substitution induced giant strain gradient and lead to the formation of dislocation. This\nprocess amplifies the substitution effect and facilitates the structural transition to the trigonal\nphase. In the CVT growth process, the system remains close to thermodynamic equilibrium,\nwith slow transport kinetics and minimal thermal fluctuations. As a result, temperature\nvariations and growth-related stress are significantly suppressed, and thus are not considered\nfactors for the dislocation formation.\nFig. R21 (In Fig. 2c in Manuscript) Geometric phase analysis (GPA) was applied to extract\nthe strain distribution at the edge dislocation regions.\nTo address this comment, we have further emphasized the gradient strain in the manuscript as\nfollows: “For Se content below 2.6%, the substituted Se atoms are randomly and sparsely\ndistributed within the 2D layers, and the induced point defects are insufficient to distort the\nlattice, maintaining a homogeneous single-phase structure. And the independent monoclinic\nphase was observed, characterized by the 98.9° angle between the a and c axis and clearly vdW\ngap (Supplementary Fig.5 and Fig. 8). With the Se content increase to x = 2.6%, local Se\nenrichment near defect sites induces sufficient strain to distort the lattice, leading to the\nformation of dislocation and trigonal phase. The coexistence of monoclinic and trigonal phase,\nalong with dislocations at the boundaries of these two phases are clearly observed in Fig. 2c.”\n6) Moreover, the authors correlated the relaxor behavior with the formation of the edge\n29\ndislocations, which again raised concern about the clear mechanism of the current work. A\nsystematic investigation was required to establish whether these dislocations were truly induced\nby substitution or were inherent to the growth process using quantitative defect density analysis.\nResponse:\nThanks for your insightful question regarding the origin of the dislocations. To reach\nsystematic investigation, all crystals were synthesized under the same CVT growth\nconditions to ensure consistency, with Se substitution being the only variable. In the\nundoped CuInP S crystal, only a uniform monoclinic phase is observed and without the\n2 6\ndislocations (Fig. R22). In contrast, the dislocations were observed once Se is introduced at\ncontent of 2.6% (corresponding to 1.4% Se content in previous manuscript). This clear\ncompositional dependence provides strong evidence that the dislocations are formatted from\nthe Se-induced strain gradients rather than from the growth process itself.\nFig. R22 Crystal structure with Se content of 0% and 2.6%. a-c, The HAADF-STEM image of\nuniform monoclinic phase at Se content of 0%: a, The atomic image viewed along the (001)\ndirection, scale bar, 1 nm. b, The atomic image viewed along the (010) direction. scale bar, 2\nnm. c, The SAED pattern along the (001) zone axis. d-f, The HAADF-STEM image of\ncoexistence of monoclinic and trigonal at Se content of 2.6%: d, The atomic image viewed\nalong the (001) direction, scale bar. 1 nm. e, The atomic image viewed along the (010) direction.\nscale bar, 5 nm. f, The SAED pattern along the (001) zone axis.\n7) The SHG results presented for the pure composition show dark regions that are not explained.\nsymmetry was broken [Adv. Mater. 35, 2210562 (2023)], the presence of dark regions in\nnominally ferroelectric domains was relatively confusing. The authors need to discuss whether\nthese dark regions arose from structural inhomogeneity, phase coexistence, or during the\nexperimental process.\nResponse:\nTo illustrate the evolution of polarization distribution with Se substitution, the SHG\nmapping intensities in the manuscript were normalized. This normalization may have led to a\nmisunderstanding that the dark regions exhibit zero SHG intensity. In fact, the quantitative\nSHG intensity analysis shows that the dark regions on the pure crystal (x=0%) still exhibit\na clear SHG signal, which weaker than that of the bright regions. After applied the unified\nintensity scale for the SHG intensity of the 0% and 2.6% (corresponding to 1.4% Se content in\nprevious manuscript) crystal, a uniformly bright SHG contrast is observed in the pure crystal,\nindicating a uniform monoclinic phase (Fig. R23a). The few scattered dark spots may arise from\nslight local structural heterogeneity or surface morphology variations introduced during sample\npreparation and measurement. For the Se content of 2.6%, the overall SHG signal is\nsignificantly reduced compared to the pure crystal, and the SHG intensity in the dark regions is\nalso non-zero (Fig. S23b). This behavior is associated with the coexisted phase structures,\nwhere the locally stabilized trigonal regions and associated dislocations possess a weaker\nnonlinear optical response due to the low symmetry structure.\nFig. R23 SHG mapping of Se content at 0% (a) and 2.6% (b) with a selected region of 10 μm\n31\n× 10 μm, scale bar: 1 μm. The plots below display the SHG intensity values corresponding to\nthe dark regions highlighted in the mapping.\nTo address this comment, we have updated the Fig. 3e by modifying the scale bar to\nenable a clear comparison of the polarization response.\n8) Peak splitting in current-electric field loops provides direct evidence of antiferroelectric\nswitching. As the switching current measurements are accessible to the authors. They should\ninclude these results to evidence the presence of antiferroelectric-like behavior by showing the\ntwo switching peak characteristics.\nResponse:\nThank you for this important comment. We agree that the current is important for analyzing\nantiferroelectric behavior. Following your comment, the current-electric field (I–E) loops was\nsupplied in Fig. R24. The typical double current peaks are clearly observed, providing direct\nevidence of antiferroelectric switching behavior. This observation is fully consistent with the\ndouble hysteresis loop in the P–E measurements, further confirming the presence of the\nantiferroelectric in the crystals with 2.6% (corresponding to 1.4% Se content in previous\nmanuscript) Se content. Moreover, as the temperature increases, the double current peaks\nprogressively merge into a single broad peak, indicating the transition from antiferroelectric to\nrelaxor ferroelectric behavior.\nFig. S24 P–E and I–E curves at phase transition from antiferroelectric to relaxor ferroelectric\nwith Se content of 2.6%.\nTo address this comment, we have updated the Supplementary Fig. 12 in the revised\n9) The P-E hysteresis loops exhibited noticeable imprint behavior, but the origin was not\ndiscussed. Given the presence of edge dislocations, it was plausible that internal fields from\ndefect dipoles or strain gradients contribute to double-like hysteresis loops, where pinched\nhysteresis loops and linked double switching current behaviors were observed with the\nformation of dislocations [ACS Appl. Electron. Mater. 6, 8328-8338 (2024)]. The authors should\nanalyze and provide further discussion on whether the imprint was correlated with\ndislocation/defect densities.\nResponse:\nWe really appreciate that your insightful viewpoint and the relevant references. We highly\nagree with your perspective that the dislocations are closely related to the imprint\nbehavior in P–E loops. According to the previously studies, (ACS Appl. Electron. Mater. 6,\n8328–8338 (2024)), the dislocations are associated with dipole pinning and strain gradients,\nwhich result in pinched P–E loops and the appearance of double switching current peaks.\nTherefore, the double hysteresis loop with imprint behavior at x = 2.6% (corresponding to 1.4%\nSe content in previous manuscript) are primarily originates from antiferroelectric phase and\ndislocation-related pinning effect. Moreover, as the Se concentration increases, the dislocation\ndensity decreases and antiferroelectric domains become further refined. As a result, the pinched\nfeature in the P–E loops gradually weaken, and the typical antiferroelectric hysteresis loop was\nobserved at Se content of 4.2% and 6.9% (corresponding to 2.5% and 4.2% Se content in\nprevious manuscript).\nTo address this comment, we add “When the Se content is below 0.8%, the CuInP (S\n2 1–\nSe ) crystals exhibit typical FE hysteresis loops, even at the low temperature of –120 °C. Once\nx x 6\nthe Se content reaches 2.6%, the pinched hysteresis loop is observed below 0 °C, indicating\nthat the substitution of Se for S atoms disrupts the uniform ordered dipoles and induces the AFE\nbehavior. This imprint behavior mainly arises from the dislocation-related pinning effect1,2.\nWith further increasing Se content to 4.2%, the more uniformly distributed Se effectively\nreleases the accumulated strain and reduces dislocation density, while simultaneously inducing\nnonpolar regions and refining local polar nanoregions. As a result, the crystal with 4.2% Se\npresents a more typical antiferroelectric-like hysteresis loop. The refinement of AFE domains\nand reduced polar phase fraction lead to decrease in polarization. As the Se content further\nincreases to 6.9%, the slim and nearly linear hysteresis P–E loops indicate that the AFE\ndomains are further refined and weaken the polarization.” in Supplementary Note 2.\n33\n10) In Raman spectroscopy (Fig. 1b), the authors claimed that the intensities of the Cu1+\nvibration mode at 315 cm-1 weaken with higher Se content, indicating that the mixed S-Se\nchalcogen environment disrupted the hopping motion of the Cu1+ cations. However, on a close\nlook, it seems like the S-P-S and Cu+1 peaks were more widened with the increasing Se\ncontent, with Cu+1 nearly comparable in intensity to the intensity of the peak at 346 cm-1.\nResponse:\nWe would like to thank you for this careful review and point out this issue. The observed\nbroadening peaks of the S-P-S and Cu1+ vibrational modes with increasing Se\nconcentration can be ascribed to the enhanced local structural heterogeneity. The partial\nSe substitution produces gradient strain, which promote the formation of dislocation and local\ntrigonal regions. These structural inhomogeneities introduce a variety of local bonding\nenvironments and enhance phonon scattering. Consequently, the Raman peaks broaden because\nthe vibrational frequencies become more dispersed in these locally disordered regions.\nWith increasing Se content, the intensities of the 315 cm–1 and 346 cm–1 Raman modes\nexhibit a clear competitive trend. The gradual decrease in the intensity of the Cu⁺ vibration at\n315 cm–1 arises from the modification and increasing disorder of the local Cu1⁺ coordination\nenvironment caused by partial Se substitution, which disrupts the Cu1⁺ hopping-related\nvibrational coherence. In contrast, the progressively enhanced peak at 346 cm–1 reflects the\nincreased formation of mixed [P (S Se ) ]4− anionic units, which exhibits stronger\n2 1−x x 6\npolarizability and produce a higher Raman intensity. These two competing effects shift the\nlattice vibrational spectral weight from the Cu sublattice toward the anionic framework. As a\nresult, the two modes reach comparable intensities, reflecting that their vibrational\ncontributions reach a balance.\nTo address this comment, we add “And the simultaneous broadening of the both 315 cm−1\nand 346 cm−1 peaks demonstrates that Se substitution enhance lattice instability, which drives\nthe formation of heterogeneous structures.” on page 3 of revised manuscript.\n11) In Figure S2, the labeling of the horizontal axis (electric field vs. voltage) was unclear.\nResponse:\nWe are grateful for your careful and professional review of the work. To clearly illustrate\nthis point, the x-axis label was revised to display temperature at 10 °C intervals as shown\nin Fig. R25. All P–E hysteresis loops were measured under an electric field of 20 kV mm–1 at\ntemperature were supplied in the supplementary Fig. 2, which was measures at 800 Hz under\nthe same electric field (20 kV mm–1).\nferroelectric behavior with temperature in the range of –120 to 25 °C. P–E loops were measured\nfrom –120 °C to 10 °C at 10 °C intervals under the same electric field of 20 kV mm–1 and\nfrequency of 100 Hz. P–E loops were measured at 25 °C under the electric field of 20 kV mm–\n1 and frequency of 800 Hz.\nTo address this comment, we have updated Supplementary Fig. 2 by display the x-axis\nto temperature and adding the room-temperature hysteresis loops with higher frequency.\n12) For fair presentation of the results, the authors should show the full XRD results in wider\nrange.\nResponse:\nIn response, the full XRD patterns from 10° to 60° have been presented in Fig. R26.\nWhen the Se content is low, the randomly substituted Se atoms induce local strain that is\ninsufficient to distort the lattice, and the diffraction peaks remains unchanged. When the Se\ncontent reaches 1.4%, the local rich Se region leads to notable lattice distortion, giving rise to\ndislocations and the emergence of local trigonal phase. The resulting local heterogeneous\nstructure causes a broadening of the (002) diffraction peak. When the Se content exceeds 2.6%\n(corresponding to 1.4% Se content in previous manuscript), the substituted Se atoms become\nmore uniformly distributed, which helps release the accumulated strain and stabilize the\ninterlayer stacking order. As a result, the crystallinity improves and the (002) peak becomes\nsharper again without a noticeable decrease in intensity. In contrast, the reflection peaks\ncorresponding to higher-order planes, such as (004), (006), and (008), gradually decrease in\nintensity with increasing Se content, which is related to local heterogeneous structure.\nfrom 0% to 6.9%.\nTo address this comment, we have updated the Fig. 2a with full XRD pattern, and we\nadd “The clear XRD diffractions peaks at 13.4°, 27.1°, 41.3° and 56.3° correspond to the (002),\n(004), (006) and (008) planes, respectively, which indicate the 2D crystal stacks along the c-\naxis (Fig. 2a).” on page 3 of revised manuscript.\nTao Yang et al synthesized layered van der Waals CuInP (S Se ) crystals with controllable\n2 1-x x 6\ncoexistence of the monoclinic and trigonal phases, which in turn leads to the appearance of\nrelaxor-like ferroelectricity in the phase coexistence region. This, indeed, can transform the\nferrielectric long-range order of CuInP (S Se ) into relaxor-type ferroelectricity with\n2 1-x x 6\nsuperparaelectric features. As evidence of relaxor-type ferroelectric properties, the hysteresis-\nless polarization switching and diffused phase transition were demonstrated. Also, HR-TEM\nconfirms the principal role of edge dislocations play in formation and stabilization of locally\ndisordered polar structures. In overall, I support authors’ claim that a phase engineering of\nlayered van der Waals ferrielectrics, opening their great potential for advanced nanoelectronics\nand related applications.\nI recommend the work for publication after a major revision. I regard that the formulation of\nthe key findings made in this work and corresponding terminology should be reconsidered and\nsignificantly softened due to the following argumentation.\nResponse:\nWe are deeply grateful for your thorough and positive evaluation of our work. Meanwhile,\nwe sincerely appreciate for your recognition on our work as “I support authors’ claim that a\nphase engineering of layered van der Waals ferrielectrics, opening their great potential for\nadvanced nanoelectronics and related applications.” According to the insightful suggestions,\nwe conducted a thorough evaluation of the core statements and made corresponding revisions\nacross the entire manuscript to ensure clarity and accuracy. We believe that these revisions have\nsubstantially improved the quality of the manuscript and that the revised version meets the high\nstandards of this journal.\n1) To the best of my understanding and knowledge CuInP (S Se ) are ferrielectrics (FI), not\n2 1-x x 6\nferroelectrics (FE), which is confirmed by extensive literature. See experimental and theoretical\ngrounds in the works: Y. M. Vysochanskii, V. A. Stephanovich, A. A. Molnar, V. B. Cajipe, and\nX. Bourdon, Raman spectroscopy study of the ferrielectric-paraelectric transition in layered\nCuInP2S6. Phys. Rev. B 58, 9119 (1998), https://doi.org/10.1103/PhysRevB.58.9119 N.\nSivadas, P. Doak, and P. Ganesh. Anharmonic stabilization of ferrielectricity in CuInP2Se6.\nPhys. Rev. Res. 4, 013094 (2022), https://doi.org/10.1103/PhysRevResearch.4.013094 A. N.\nMorozovska, E. A. Eliseev, A. Ghosh, M. E. Yelisieiev, Y. M. Vysochanskii, and S. V. Kalinin.\nAnomalous Polarization Reversal in Strained Thin Films of CuInP2S6, Phys. Rev. B 108,\n054107 (2023) https://link.aps.org/doi/10.1103/PhysRevB.108.054107. The AFE-FI-FE\ntransitions of hysteresis loops are very clearly seen from Supplementary Fig. 2, that in my\nopinion, can be moved to the main text. Also, the figure shows a very significant contribution\nfrom the electric conductivity (losses) to the hysteresis loops.\nWe appreciate for point out this important point, which is really helpful to improve the\naccuracy and significance of our conclusion. According to our analysis, the CuInP 2 (S 1-x Se x ) 6\npossesses a ferrielectric (FiE) configuration, arising from opposite displacements of Cu\nand In ions with unequal magnitudes, leading to a nonzero net polarization. During the\npolarization switching process, the reversal is predominantly governed by the displacement of\nCu ion, while the contribution of In ion is negligible. As a result, it macroscopically exhibits\nferroelectric (FE-like) behavior characterized and a typical P–E hysteresis loop.\nAs you mentioned, the intermediate state may appear during the transition from AFE to\nFE phase. While this state with pinched P–E loop exhibits the profiles typical of imprint\nbehavior. This behavior mainly arises from the antiferroelectric phase and dislocation-related\npinning effect as shown in P–E curves in Supplementary Fig. 2 and has been discussed in our\nresponse to Reviewer 3, Question 9. The exam of the I–E loops at low temperature for 2.6%\n(corresponding to 1.4% Se content in previous manuscript) Se exhibits a typical double-peak\nfeature, which is inconsistent with the three current peaks for ferrielectric behavior with triple\nhysteresis loop (Fig. R27). In addition, to eliminate the influence of leakage current, we have\nsupplemented the P–E loops measured at 800 Hz at room temperature in Supplementary Fig. 2.\nFig. R27. P–E and I–E curves at phase transition from antiferroelectric to relaxor ferroelectric\nwith Se content of 2.6%.\nTo address this comment, we have provided additional clarification in revised manuscript\nas “Notably, the crystal of CuInP (S Se ) possesses a ferrielectric configuration, which the\n2 1-x x 6\nCu and In ions exhibit antiparallel displacements with unequal magnitudes, resulting in a non-\nzero net polarization33. Since the contribution of In ion to the polarization is relatively weak,\n38\nthe ferroelectric behavior is fundamentally governed by the displacement and positioning of Cu\nions,” on page 5 of revised manuscript. We have updated the Supplementary Fig. 2 with P–E\nloops measured at 800 Hz at room temperature to minimize the influence of leakage current.\n2) The ferrielectric 3-well or 4-well free energy potential can be transformed into the relaxor-\nlike relief, but the relative dielectric permittivity (shown in Fig. 3e) is not very high (or “giant”)\nas it should be in true relaxors. The diffuse maximum of permittivity does not exceed the values\nof 200-350, being much smaller than 104-106 (as in ferroelectric relaxors). While the frequency\ndispersion of the permittivity is present in the frequency range 300 Hz – 3 MHz, the dispersion\nis much weaker than in ferroelectric relaxors. As a matter of fact, the ferrielectric with a diffuse\nphase transition and relaxor-like properties is obtained.\nResponse:\nWe really appreciate your valuable comment, which provides an accurate description for\nthe relaxor-like behavior. As you noted, the relatively small frequency dispersion indicates that\nthe relaxor-like ferroelectric is observed, and the corresponding revisions have been made in\nthe manuscript.\nTo address this comment, we add “In contrast to the sharp dielectric peaks for x = 0%\nwith typical ferroelectric phase, the gradual shift of phase transition temperature at relatively\nsmall range with frequency confirms the relaxor-like ferroelectric at composition x = 2.6% and\n6.9% (Fig. 3e).” on page 6 of revised manuscript.\n3) I failed to understand why the author called the behavior shown in the bottom of Fig.3e as\nsuperparaelectric (SPE). Here the relative dielectric permittivity is relatively low; it increases\nabove 40-50 °C. The increase of dielectric permittivity above 40-50 °C can be related to the\nthermally-activated conductivity (Maxwell-Wagner type effects). This is clearly seen from\nSupplementary Fig. 8 and Supplementary Fig. 9, where the imaginary part of dielectric\npermittivity and tangent of dielectric losses increases strongly above 40-50 °C. True relaxors\nand superparaelectrics (if any exist) should be dielectrics in the actual range of these properties’\nobservation. Also, the polarization in saturation should be very high in a relaxor that response\nto the external electric field as a super-dipole with a large dipole moment. However, the\npolarization shown in Supplementary Fig. 2 for x=4.2 % is relatively low (below 2 µC/cm2),\nThus, authors rather observe some SPE-like features, but not the SPE state.\nResponse:\nThank you for your careful analysis and we fully agree that the SPE-like state is observed\nrather than the SPE state. As the Se content increases, the polarization gradually decreases and\nreaches below 2 µC cm–2 at 6.9% (corresponding to 4.2% Se content in previous manuscript).\nThis reduction mainly originates from the decreased fraction and refinement of the local polar\nregions, as evidenced by the PFM images. The relatively low polarization also indicates that\nresidual polar regions still remain and the system has not completely transformed into a\nparaelectric phase. Therefore, we define this intermediate transition state as SPE-like state\nin this manuscript.\nTo address this comment, we add “which displays slim and nearly linear hysteresis loops\ncharacteristic of the SPE-like state.” on page 5 of revised manuscript.\n4) The work is very well illustrated, but at the same time, the important detail shown in the\nfigures are not commented neither in the main text nor in the Supplement.\nResponse:\nThanks for your positive evaluation of our illustrations. In response, we have carefully\nrevised the manuscript to ensure that the details highlighted in the figures are clearly described.\nSpecifically, we supplied additional analyses and explanations for Fig. 1b, Fig.2a,\nSupplementary Fig. 2, Fig. 3e and Supplementary Fig. 3 in the main text or in the\nSupplementary Information. These improvements enhance the clarity and enable a more\ncomprehensive understanding of the data presented.\nTo address this comment, we add analysis of Raman spectroscopy in Fig. 1b as “And the\nsimultaneous broadening of the both 315 cm−1 and 346 cm−1 peaks demonstrates that Se\nsubstitution enhance lattice instability, which drives the formation of heterogeneous structures.”\non page 3 of revised manuscript.\nWe add analysis of XRD in Fig.2a as “As the Se content surpass 2.6%, the more uniform\ndistribution of Se effectively releasing the accumulated strain and improving the local stacking\norder, resulting in the recovery of sharper (002) diffraction peak. And the diffraction peaks shift\nto larger angles, reflecting a reduction in interlayer distance along the c-axis. Additionally, the\ngradual decrease in intensity of the (004), (006) and (008) peaks with increasing Se content is\nrelated to local heterogeneous structure.” on page 4 of revised manuscript.\nIn addition, we add analysis of Supplementary Fig. 2 in Supplementary Note 2 as “The P–\nE loops from –120 °C to 10 °C were recorded at intervals of 10 °C under an electric field of 20\nkV mm–1 and frequency of 100 Hz. The room temperature P–E loops were measured under the\nsame electric field but at a higher frequency of 800 Hz to minimize the influence of leakage\ncurrent. When the Se content is below 0.8%, the CuInP (S Se ) crystals exhibit typical FE\n2 1–x x 6\nhysteresis loops, even at the low temperature of –120 °C. Once the Se content reaches 2.6%,\nthe pinched hysteresis loop is observed below 0 °C, indicating that the substitution of Se for S\natoms disrupts the uniform ordered dipoles and induces the AFE behavior. This imprint\nbehavior mainly arises from the antiferroelectric phase and dislocation-related pinning effect1,2.\nWith further increasing Se content to 4.2%, the more uniformly distributed Se effectively\nreleases the accumulated strain and reduces dislocation density, while simultaneously inducing\nnonpolar regions and refining local polar nanoregions. As a result, the crystal with 4.2% Se\npresents a more typical antiferroelectric-like hysteresis loop. The refinement of AFE domains\nand reduced fraction lead to decrease in polarization. As the Se content further increases to\n6.9%, the slim and nearly linear hysteresis P–E loops indicate that the AFE domains are further\nrefined and weaken the polarization. In addition, for Se content of 2.6% and 4.2%, the\nantiferroelectricity gradually transforms into RFE as the temperature increase, which is\nattributed to the thermal activation causes the antiferroelectric domains convert into local\nordered polar nanoregions. These polar nanoregions create numerous defects and domain walls,\nleading to a significant increase in leakage current. Therefore, the P–E hysteresis at room\ntemperature loops were measured at a higher frequency (800 Hz) to effectively minimize the\nleakage current effect. With increasing Se content, the gradually weakened hysteresis and\nreduced maximum polarization at room temperature confirms the transition from FE to RFE\nand SPE behavior.”\nWe add analysis of Fig. 3e and Supplementary Fig. 3 in Supplementary Note 3 as “In the\nSe-doped CuInP (S Se ) crystals, the increased electrical conductivity caused by Se\n2 1–x x 6\nsubstitution introduces significant Maxwell–Wagner polarization at low frequencies. At the low\nfrequency, the charge carriers can accumulate at internal interfaces or defects, producing an\nadditional space-charge polarizations. In this case, the overall dielectric response arises from\nthe superposition of dipolar and space-charge polarizations, and the relaxation process\ncorresponding to the maximum dielectric constant occurs at relatively higher temperatures. As\nfrequency increases, this space-charge contribution is gradually suppressed because carriers\ncannot follow the rapid electric field. The dielectric response then becomes dominated by faster\ndipolar relaxation, which maximum dielectric constant occurs at lower temperatures. Thus, the\ndielectric peak shifts toward lower temperature with increasing frequency.”\nIn summary, we would like to express our sincere gratitude for your time and effort in\nthoroughly considering our manuscript. These insightful comments and constructive\nsuggestions help us further improve our work's quality and clarity. We believe the revised\nmanuscript can meet the high standards of this authoritative journal. Thank you very much!\n(Manuscript No: NCOMMS-25-70669B)\nReviewer #1 (Remarks to the Author):\nThe authors have made substantial improvements to the quality of the work by adding\nmore details and explanations, especially related to the growth and structural characterizations,\nwhich were lacking in the original manuscript. I appreciate the authors' detailed responses and\nconducting further experiments to clarify the issues raised by the reviewers. I recommend the\narticle for publication.\nHowever, I do want to point out one misunderstanding in one of my questions in the first\nround of reviews related to the measurement of the strain in relaxors. PFM, especially when\noperating near the resonant enhanced mode, is not an accurate method to quantify the strain as\nthe PFM amplitude signal can be convoluted by many different factors (see Neumayer et al.,J.\nAppl. Phys. 128, 171105 (2020)). For an artifact free electromechanical strain measurement,\ntypically laser Doppler vibrometer setups are used to generate butterfly-type strain curves (see\nFig 1 in Pan et al., Nat Mat, 23, 944-950, 2024) which can be used to accurately quantify the\nstrain change. However, this point does not affect the main message of the manuscript.\nResponse:\nWe would like to express our sincere appreciation for your positive evaluation of our\nrevised manuscript and for recommending our article for publication. Thank you for pointing\nout that strain quantification using PFM, especially under resonance enhanced mode, is not\naccurate. We also learn that laser Doppler vibrometer (LDV) measurements can provide a more\nreliable and quantitative approach to obtain butterfly-type strain curves. This valuable\nsuggestion will greatly benefit our future studies, and we sincerely appreciate your guidance.\nReviewer #2 (Remarks to the Author):\nAfter carefully reviewing the responses by the authors for the given comments, I find that\nthey have provided satisfactory explanations and adequately addressed all the concerns raised.\nThe revisions have significantly improved the clarity and quality of the manuscript. Therefore,\nI recommend the manuscript for acceptance in the journal.\nResponse:\nThank you for your time and effort in carefully reviewing our manuscript, we greatly\nappreciate your constructive feedback and your recommendation for acceptance.\nThe Authors have revised the manuscript and thus improved the quality of the manuscript,\nthus, it is considered for publication to Nature Communications.\nResponse:\nWe really appreciate for your positive evaluation of our revised manuscript. Thank you for\nyour recommendation for acceptance.\nReviewer #4 (Remarks to the Author):\nAuthors answered to all my critical remarks and suggestions in a constructive way. I\nrecommend the revised work for publication.\nResponse:\nThank you for your careful review. We appreciate your recommendation for publication.\n2"}]

Decision letter