1. INTRODUCTION

Ultrafast laser micromachining has emerged as a versatile and high-precision platform for three-dimensional (3D) micro- and nanofabrication in transparent materials [1]. Enabled by nonlinear absorption phenomena, such as multiphoton absorption and avalanche ionization [2,3], femtosecond laser direct writing (FLDW) allows the fabrication of complex internal structures with submicron resolution, without the need for masks or chemical processing [4–7]. This capability has fueled advances in a range of applications, including optical waveguides writing [8–15], stealth machining [16,17], ferroelectric domain engineering [18–20], multiphoton polymerization printing [21–24], and functional materials engineering [6,7]. Apart from this, another key strength of FLDW lies in its ability to fabricate arbitrarily designed 3D structures in a single step without intermediate alignment or layer-by-layer stitching. This feature is particularly critical for emerging applications, such as high-dimensional topological photonic chips [25,26], fiber-to-photonic integrated circuits (PICs), and PIC-to-PIC high-dimensional optical interconnects [27–29], where precise, continuous, and spatially variant geometries are essential. In these systems, the ability to directly inscribe complex waveguide trajectories, including helical, curved, or multidepth crossings, provides unprecedented design freedom and integration flexibility.

Nevertheless, fabricating continuous 3D structures with high uniformity and consistency over large depth changes [7] remains a persistent challenge. As the writing depth increases, the interaction between the laser and matter becomes progressively distorted due to spherical aberration, scattering, and nonlinear energy attenuation [30], particularly when the laser beam traverses through material interfaces [e.g., air–glass boundaries [31–33], Fig.  1 (a)], inhomogeneous media [e.g., multiphoton 3D printing photoresists [34], Fig.  1 (b)], or heterogeneous and multilayer substrates [e.g., anisotropic 2D materials [35] or monolayer stacks [36], Fig.  1 (c)]. These distortions impact both the optical focus profile and the thermal accumulation behaviors, leading to depth-dependent variations in feature dimensions, refractive index contrast, and optical loss, all of which degrade the performance and yield of photonic devices [Fig.  1 (d)].

Fig. 1. (a)–(c) are the refractive index interface present in the laser path during micromachining. (d) and (e) are the schematic of the manufactured result using different speed control method. (f) The side view of the waveguide that dived from 41 to 507 µm beneath the substrate surface with the white dashed line indicating the surface position.

More critically, nonlinear absorption further amplifies these distortions: since multiphoton absorption requires extremely high intensities, even small fluctuations in energy can lead to substantial variations in the effective focal intensity, resulting in large discrepancies in the material response. Several strategies have been proposed to mitigate focal distortions at specific depths, including adaptive phase correction [37–40], in situ depth detection and correction [41], and defect-compensated methods [42]. However, these approaches are generally limited to discrete depth compensation and lack the flexibility to support the fabrication of continuous 3D line-like structures. Maintaining uniformity and consistency in such structures, especially when they extend across large depth ranges, remains a fundamental challenge in deep 3D micromachining. Therefore, addressing depth-dependent energy inconsistencies is vital to unlocking the full potential of FLDW for next-generation, large-scale photonic systems.

Under thermal FLDW conditions [43,44], in addition to energy field distribution in the focus, it is of paramount importance to consider the absorbed energy per unit exposed volume in the context of thermal accumulation in femtosecond laser fabrication, where the thermodynamic effect plays a pivotal role in influencing the consistency of the machining process [45,46]. In this case, the net energy flux, accounting for both the single-pulse energy distribution and the thermal accumulation effect induced by multiple pulses, can be defined as follows:

Fig. 2. Schematic of the machining system and the technical flowchart of the speed control method. With the assistance of the waveguides jiggle data extracted by digital image processing and the position error data from the stage A3200 Digital Scope software, the positional error introduced was evaluated for optimization.

where the  is the beam waist diameter, the  is the average single pulse flux affected by the energy distortions, the  is the repetition, and the  is the scan speed [47]. When using high-repetition-rate pulsed lasers satisfying the condition:

the laser pulse can be treated as a continuous heat source [48]. In this regime, the thermal accumulation process is governed by heat conduction dynamics in solids [49], and  plays a critical role in modulating this effect. Unlike previous efforts, including adaptive phase correction [37–40], energy correction [42] and depth correction [41], which primarily corrected , our approach leverages  to compensate for variations in the absorbed energy flux. Based on this principle, we developed a novel adaptive energy-flux compensation (AEFC) method for continuous structure processing [Fig.  1 (e)]. The method synergistically controls energy depositions and heat accumulation in real time by dynamically adjusting the motion stage precisely at a high frequency of 1–2 kHz without introducing additional optical path settings. This represents an order of magnitude higher than conventional adaptive phase correction based on spatial light modulators. In principle, the consistent structural processing of materials with variations in depth of up to µ is a feasible undertaking. Comparing to other depth compensation method, the AEFC method is a universal adaptive compensation with high modulation, which can function in combination with other beam-shaping method (see Table S1 of Supplement 1, Appendix A). In experimental trials, the slit beam-shaping method in combination with the AEFC method was utilized to fabricate waveguides achieving an insertion loss () that remained consistently over a depth variation of 500 µm [Fig.  1 (f)].

2. RESULTS

A. Setup and Measurement

As described in Eq. ( 1 ), the scan speed varies along the machining route that decides the energy distortions. However, the variation of the  with depth is critically complex, so that the curve function of the scan speed, which is controlled for energy-flux compensation, is not simple straight line or sine curve. Consequently, the ramp control function of the motion stage is unable to fulfill the requirement of the complex scan speed curve. In such cases, the energy flux, which varies with depth, necessitates compensation by a complex speed curve that must be regulated at the microinterval level for the continuity of structures. To solve this problem, we proposed specifying the scan speed in each microinterval to fit the complex scan speed curve during the machining. In order to provide a concrete illustration of this speed control method, consideration may be given to 3D waveguides that diving through large depth change. As illustrated in Fig.  2 , the sample is moved with the motion stage during the machining process. Concurrently, the focus of the femtosecond laser undergoes relative movement within the sample. The positional deviation of the motion stage is collected from the feedback of the controller using the A3200 Digital Scope software. In parallel, skeleton extraction of image processing technique [50] is employed to extract the skeleton of the waveguide (more details in Supplement 1, Appendix D). This enables the geometry of the fabricated waveguide to be quantized from the microscope image. The data collected from both sources are utilized for the purpose of evaluating and optimizing the speed control method.

Fig. 3. Captured waveguide skeleton (blue) and the stage Y-axis position error (red) at the corresponding moment.

To verify the feasibility of controlling the speed in each microinterval, we fabricated straight waveguides with an increasing scan speed from 8 to 80 mm/s, which has a length of 9 mm, and the interval is fixed as 0.1 mm. The waveguides are jiggled as shown in the microscopic top view of Fig.  2 . In order to find the cause of the problem, we analyzed the skeleton of the 0.5–1 mm portion (a low scan speed corresponding from 8 to 12 mm/s) of waveguides manufactured at different moments, surprisingly found that they match up well (as shown in Fig.  3 ). Meanwhile, the Y-axis position error curves of the stage, which correspond to the waveguides, match up well too, and the spacing of peak to peak corresponds to the control interval 0.1 mm. It indicates that the intense jiggle is not random but strongly relate to the control of the stage. As such, we can optimize the control method to achieve the scan speed as we wish while eliminating the jiggle.

Fig. 4. (a) Schematic of the nonadaptive control method of scan speed. (b)–(d) The scatter diagrams of the motion stage position error at various fixed  and speed zones. (e) The schematic of the adaptive control method. (f) The scatter diagrams of the motion stage position error with adaptive  and various speed zones.

B. Speed-Adaptive Sampling Method

As illustrated in Fig.  4 (a), the microinterval is labeled as , which is fixed as a constant while manufacture processing. The first row of Fig.  4 (b) shows the scatter diagram of the position error of X axis and Y axis during the machining processing with the stage accelerating from 16 to 24 mm/s at 1 mm length, and the  is 0.1 mm. It is evident that, although the mean of the error is correct, the variance is large. This will induce jiggle with the waveguides. The alteration in the speed-up zone, from 16–24 mm/s, 40–48 mm/s to 64–72 mm/s, resulted in a reduction in the variance of the position error. The results demonstrate that the Y-axis position error of stage diminishes as the speed increases, while the  is 0.1 mm. When we reduce the  from 0.1 mm to 0.05 mm, the variance is diminished overall, but the variance of the low scan speed portion (16–24 mm/s) remains unacceptably high [as shown in Fig.  4 (c)]. In order to investigate the limits of , it was set at 0.01 mm, a value that represents a significant reduction from previous. As illustrated in Fig.  4 (d), the predicted reduction in the variance of the low scan speed portion was observed, and the mean was found to be correct. Conversely, the variance of X-axis position error at the high scan speed portion (40–48 mm/s and 64–72 mm/s) increased significantly, indicating that the standstill of the X-axis movement occurred during the actual processing, which resulted in discontinuities in scanning speed. This, in turn, affected the uniformity of thermal accumulation, thereby creating scattering points and dimensional inhomogeneities in the waveguide (see Fig.  2 ). These findings indicate that distinct  are applicable to varying speeds. A simple way to optimize the  is to make it vary with the scan speed so that the match between the two can be set in the relation mentioned above [as shown in Fig.  4 (e)]. So, we transfer the fixed control interval in space to time (more details in Supplement 1, Appendix E). Then the actual  is expressed as

where  is the scan speed at the position of interval beginning. This evolved , designated the adaptive interval method, is characterized by a variable speed-dependent nature and is deemed suitable for implementation with any specific scan speed curve varying with the position of X. In contradistinction to the fixed interval method, the adaptive control method, when the  was set at 0.5–1 ms, was found to be capable of maintaining reduced variances and tolerable means at a variety of speed zones [see Fig.  4 (e)]. This finding indicates that the adaptive interval method can not only reduce the position error of the Y axis, but also eliminate the standstill of the X axis.

To validate whether the adaptive control method could mitigate jiggle at low scan speeds and avoid scattering points at high speeds, we fabricated a series of waveguides. The scan speed was set from 8 to 80 mm/s, and the X range was 9 mm. As previously described, the skeletons of the waveguides were obtained through the use of skeleton extraction of image processing techniques. In order to assess the jiggle of the waveguides in comparison to the ideal straight waveguide, it was necessary to consider the tilt of the sample during the taking of micrographs. A linear fit was initially generated for the skeleton, and then the skeleton was compared to the fitted line, calculated variance () of the deviation of each pixel from the fitted line. The  is described by the formula:

where  is the deviation of each pixel from the fitted line, and  is the numbers of the pixels. The heatmap of  that varies with the scan speed and control interval is shown in Fig.  5 (a). The red boxed area means that the matching of the speed and the control interval will generate the scattering point although the jiggle of the waveguide is small enough. There is a contiguous area [the blue dashed line boxed area in Fig.  5 (a)] without the intense jiggle and scattering points in the heatmap. Consequently, it can fabricate a nearly straight waveguide with an arbitrarily variable speed when the  was set as 0.5–1 ms. Meanwhile, the insertion loss of the tapered waveguides, fabricated with adaptive interval and various fixed interval, was measured. As the fixed interval is increased, the insertion loss initially decreases and then increases due to the gradual reduction in the scattering point and the subsequent gradual increase in the tremors of the waveguides. The waveguide with adaptive interval exhibits a significantly reduced insertion loss (1.07 dB) in comparison to the fixed.

Fig. 5. (a) Calculated  of the waveguide indicating the positional deviation from the ideal straight line. (b) The insertion loss of the 1 cm length sample manufactured by using various control intervals.

C. Energy-Flux Compensation Principle

Up to this point, we have focused our efforts on optimizing the control of the scan speed and improving the precision of micromachining. Based on our method, we have also sought to exert control over the energy flux, which affects the thermal accumulation. The FLDW technology provides a classic example of the application of thermal accumulation. In order to enhance the density of the waveguides array, it is necessary to utilize multilayer waveguides, which include the waveguide coupling between the different layers and cross several layers (namely large depth), which are designed for specific functions. These structures mean that the depth varies continuously. However, the energy intensity at the objective focus is found to vary with depth change, as a consequence of the spherical aberration generated by the air–material interface. The spherical aberration is expressed using the analytic function [32]:

where  and  are the refractive indices of the air and materials, respectively, and  is the normalized objective lens pupil radius. Considering the multilayer materials, the spherical aberration is only subject to variation with the depth and the mismatch of refractive index between the air and the fabricated layer. It has been established that other layers induce a fixed spherical aberration, and that they can be regarded as parallel plates (more details in Supplement 1, Appendix F). It can thus be concluded that the spherical aberration induced by multilayer materials (namely, multirefractive interfaces) can be simplified into a single interface. Subsequently, the three-dimensional Fourier transform [51] was used to calculate the field distribution in the focal region affected by the refractive index mismatch interface. The energy deposition is influenced by the intensity distribution, which in turn affects the thermal accumulation. This is considered the key factor in determining the diameter of the waveguides. In order to achieve the consistent performance of the waveguides that links several different layers, it is necessary to compensate the effect of the spherical aberration on the energy flux so that the waveguides diameter can be controlled in the same size. Here, we demonstrate the influence of spherical aberration on a  microscope objective (Olympus UPLFLN, ) that can be applied for various coverslip thicknesses from 0.11 to 0.23 mm, and verify the feasibility of compensating for the energy-flux change caused by spherical aberration through controlling the scan speed within a certain depth range.

When the coverslip thickness was set at 140 µm, which effectively corrects the spherical aberration at this depth, the simulation of the field distribution with the Y-Z plane after the astigmatic laser generated by 0.4 mm slit was focused into the materials by the ideal microscope objective (, ) at various depths is shown in Fig.  6 (a). It is evident that the energy intensity at the corrected depth is the greatest, with the others gradually decreasing as the depth is further from 140 µm. This is due to the fact that the aberration increases with depth change. Concurrently, the area that has been normalized with respect to the energy threshold and the included energy flux were calculated, as illustrated in Fig.  6 (b). Both of these quantities exhibit a similar trend, whereby they decrease as the distance from the spherical aberration-corrected depth increases. This indicates that the diameter of the waveguides will also decrease in a similar manner, which is fabricated based on the thermal accumulation. Concurrently, experiments were conducted to compare the simulation. The incident astigmatic laser was focused by the  microscope objective, with a coverslip thickness of 140 µm. The cross-sectional micrograph of waveguides at varying depths, fabricated with the same incident pulse energy 0.92 µJ and scan speed 20 mm/s, is presented in Fig.  6 (c). This figure clearly demonstrates a correlation with the simulation presented in Fig.  6 (a), with the largest waveguide diameter occurring at the corrected depth. Furthermore, the waveguide matrix was fabricated with varying depth, scan speed, and the same incident pulse energy 0.92 µJ, and the diameter drawn as a contour map was measured [shown in Fig.  6 (d)]. It is evident that the waveguides diameter influenced by the spherical aberration can be compensated by control the scan speed fitted to the contour in Fig.  6 (b) from depth 60 to 220 µm. In order to achieve larger compensable depth, we test another industrial objective (Nikon LU Plan ELWD , , ) with a smaller NA and longer work distance, shown as Fig. S6a, from depth 50 to 540 µm. The observed change trend in the diameter of waveguides is somewhat peculiar. The industrial microscope objective does not correct for the spherical aberration associated with the thickness of the coverslip, and the diameter should decrease monotonically as the depth increases. However, the diameter demonstrates a maximum value at a depth of 230 µm. This is because the slit was positioned at a limited distance from the pupil, which introduced negative spherical aberration to the infinity-corrected objective lens [52]. However, the positive spherical aberration induced by the air–material interface serves to compensate for this at a certain depth (more details in Supplement 1, Appendix G).

Fig. 6. Results of the simulation and experiment investigating the influence of spherical aberration on a  microscope objective. (a) The Y-Z energy distribution simulation of astigmatic beam at various depths. (b) The energy flux varying with the depth. (c) The cross section of waveguides fabricated at various depths shows the impact of the energy flux. (d) The contour map of waveguides diameter varying with depth and scan speed.

3. APPLICATION: DEEP FLDW FOR MULTILAYER 3D PHOTONICS CIRCUITS

The AEFC approach, when combined solely with either high-order vector optical fields [28,53] or multifoci-shaped beams [54], enables one-pass fabrication of continuous structures. This combination maintains consistent performance over depth variations exceeding 500 µm, achieved exclusively by compensating for the positive spherical aberration at the air–material interface. Furthermore, we propose that the aberration correction plane (or zero-aberration plane, ) can be positioned deep within the material using correction mechanisms such as SLM-based wavefront correction [40,55] or objectives with correction rings [56,57] [Fig.  7 (a)]. This enables AEFC to dynamically compensate for both positive and negative spherical aberrations, allowing for consistent fabrication across a total depth variation exceeding 1000 µm (i.e., µ). To demonstrate the performance of the AEFC, we fabricated the waveguides by employing the fixed scan speed and AEFC separately, with the spherical aberration correction introduced by the slit beam-shaping method. The spherical aberration corrected plane was set about 230 µm under the surface. The waveguide ports were set at 40–510 µm, which correspond to deep changes in the platform from 10 to 310 µm increasing with equal intervals, as shown in Fig.  7 (b). The astigmatic laser generated by the 0.4 mm slit was focused in the commercial glass Eagle XG with a length of 2 cm by the  microscope objective with the stage moved at the scan speed achieved from Supplement 1, Appendix H. Compared to the fixed speed method (the waveguide cross section shown in Supplement 1, Appendix I), the AEFC method can maintain the waveguide diameter regardless of different depth with the assistant of energy-flux compensation, which is evident in Fig.  7 (c). Thus, the output modes’ profile of the waveguide with different depth changes is similar [as shown in Fig.  7 (d)]. Subsequently, the insertion loss at 1550 nm of the waveguides with varying depth was evaluated. In the four groups of samples that we have tested, the average insertion loss at each depth was 0.878 dB at maximum and 0.744 dB at minimum. The average insertion loss value for all depth was 0.803 dB. Therefore, the variation in insertion loss is estimated in the range of . Meanwhile, the propagation loss is estimated 0.26 dB/cm @1550 nm (measurement details in Supplement 1, Appendix J). Based on previous work [58], the refractive index contrast of waveguides written by slit beam-shaping method in Eagle XG is .

Fig. 7. (a) Principle of the combination the aberration correction and AEFC for machining whose depth change over 1 mm. (b) The schematic of the structure. (c) The cross section of waveguides fabricated with the AEFC increasing by the same increment depth. (d) The output mode profile of the waveguides with the fixed speed and the AEFC. The scale bar is 10 µm. (e) The insertion loss at 1550 nm of the waveguides with different depth change.

4. CONCLUSION

This study addresses the challenges of focal field distortion and energy attenuation encountered during deep three-dimensional structuring using ultrafast laser processing in transparent materials. An adaptive energy-flux compensation (AEFC) method is proposed, which establishes a depth-dependent mapping between optical field distortion and material response. By implementing a coordinated regulation strategy between scanning speed and energy input, the method enables precise, dynamic control of spatiotemporal energy distribution. Experimental results demonstrate that the AEFC approach achieves real-time depth compensation exceeding µ in continuous structures, with a modulation rate increased to 1000–2000 Hz, significantly enhancing both processing efficiency and uniformity. The proposed method overcomes the limitations of traditional adaptive compensation techniques in scenarios involving rapid depth variation, providing a novel solution for high-precision 3D laser fabrication in transparent materials. The spatiotemporal energy co-regulation mechanism demonstrated by AEFC offers strong generality and scalability, establishing a key technological foundation for the development of large-scale integrated 3D photonic chips and glass-core optoelectronic co-packaging platforms. This work holds substantial promise for both engineering applications and further research exploration.